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PROTECTION OF THE SOLAR POWER SATELLITE 


Abstract 


This report examines theoretically several features of the interactions 
of the Solar Power Satellite (SPS) with its space environment. We calculate 
the voltages produced at various surfaces due to space plasmas and the 
plasma leakage currents through the kapton and sapphire solar cell blankets. 
At geosynchronous orbit (GEO), this parasitic power loss is only 0.7$, and 
is easily compensated by oversizing. At low-earth orbit, (LEO), the power 
loss is potentially much larger (3$), and anomalous arcing is expected for 
the EOTV high voltage negative surfaces. Preliminary results of a three- 
dimensional self-consistent plasma and electric field computer program are 
presented, confirming the validity of the predictions made from the one- 
dimensional models. Lastly, the report considers magnetic shielding of the 
satellite, to reduce the power drain and to protect the solar cells from 
energetic electron and plasma ion bombardment. We conclude that minor modi- 
fications can allow the SPS to operate safely and efficiently in its space 
environment. The SPS design employed in this study is the Jan 25, 1978 MSFC 
baseline design utilizing GaAs solar cells at CR-2 and an aluminuum structure. 
Subsequent design changes will substantially alter the basic conclusions in 
this report. 



Introduction 


Space Is by no means empty. It contains light, magne- 
tic fields and both neutral and charged particles. The light 
energy is the raison d'etre for space power generation; but 
it can also eject photoelectrons from satellite surfaces, 
giving the surface a positive charge and giving it an effec- 
tive conductivity (Pelizzari and Criswell, 1978). 

Magnetic field strene*‘is in the earth's vicinity range 
from 6 x 10“ S T (0.6) Gauss) at the earth's poles to 2 x 
10 _T ( 2y ) in the neutral sheet in the magnetotail ( lY 
o 10” s Gauss). At the geosynchronous orbit, the magnetic 
field strength is roughly 1 x 10“ 7 T (100 Y ). A magnetic 
field of this strength causes no threat per se to spacecraft 
operations; however, it plays a fundamental role in trapping 
energetic particles. These trapped particles respond not 
only to the Earth's magnetic field, but spacecraft fields as 
well, especially for spacecraft large in comparison to par- 
ticle gyroradii (Reiff, 1976; Reiff and Burke, 1976). 

Neutral particles have little effect rn spacecraft 
operations above ^ 600 km; however, neu- 
trals can charge-exchange in the EOTV thruster beam (see 
below) . 

Charged particle populations at synchronous orbit are 
jf several types and are illustrated in Figure 1. The inner- 
imost region is the plasma sphere, a torus-shaped locus of 
relatively dense ( ^ 100/cm 3 ), cool (kT 1 eV) plasma that 
has evaporated from the ionosphere. Because of the low 
energies of the pi asmaspheric ions, they are considered harm- 
less (Reasoner et al . , 1 976 ) ; however, they can be acceler- 
ated by spacecraft electric fields to energies high enough 
to do damage (tens of kilovolts). Imbedded in the plasma- 
sphere are the radiation belts, regions of very low density 
but quite high energy (tens to hundreds of kilovolts) trap- 
ped radiation. This radiation can cause hazards to men and 
solar cells. 

The remaining plasma population that can penetrate to 
geosynchronous orbit is the plasma sheet (Fig. 1). This 
tenuous plasma (0.1-1/cm 3 ) is considerably warmer (kT on the 
order of kiloelectron volts) than the plasmasphere (Garrett 
and DeForest, 1979). In addition, its presence at geosyn- 
chronous orbit is *,ssociated with substorm activity, when 
both the fluxes and energies are higher. It is this kind of 
plasma that contributes most strongly to spacecraft charging 
and its concomitant disruption of satellite systems (Inouye, 
1976). 



This report concentrates on spacecraft charging and its 
effects on solar power satellite (SPS) systems, in parti- 
cular the NASA/Marshall Space Flight Center (MSFC) baseline 
design (Hanley, 1978). "Worst case" plasma environments are 
used to determine possible charging hazards. Spacecraft 
charging is the principal focus of this paper since its 
effects can he severe: arc generation from exceeding break- 

down voltages, direct electrical component damage from tran- 
sients, disruption of logic and switching circuits from 
electromagnetic interference, change of reflective or ther- 
mal control surfaces due to the attraction of outgassed con- 
taminants or pitting, and shock hazards for extravehicular 
and docking activities (see DeForest, 1972; Pike and Bunn, 
1976; Shaw et al., 1976). 

We will show that, under substorm conditions, the kapton 
substrate contemplated for use as a support blanket for the 
reflectors and solar cells will be subjected to near-break- 
down voltage. Additional kapton insulation seems unfeasible 
because of weight considerations. The alternatives, higher 
conductivity substrates or conducting leads to the surfaces, 
seem more reasonable since the resulting parasitic currents 
are not excessive. The paper also will discuss the optimum 
point for grounding the spacecraft to the solar panels and 
outlines a method of using judicious routing of bus-bar 
currents to shield the satellite from particle bombardment. 
Although it is possible to use a similar method to magneti- 
cally align the satellite with the Earth's magnetic field 
(counteracting gravity-gradient torques), the fields required 
ieem unreasonably large. 

Spacecraft Charging 

A body immersed in a plasma will acquire a net charge 
from unequal fluxes of plasma particles. For most plasmas, 
the electron and ion densities N e and Nj are roughly equal, 
and the electron and ion temperatures Te and Tj are co mpara - 
ble. Thus the electron flux Je (proportional to N e /kT e /f^) 
is generally much larger than the ion flux Jj , and the body 
acquires a negative charge sufficient to bring the currents 
into balance. For stationary, isothermal, singly-charged 
plasmas, f he equilibrium unlit body potential is roughly 
(kTe/e)ln(J e /aJi ) (Whipple, 1965), where a is a parameter 
(of order unity) depending on the thickness of the sheath. 
Exposing the body to sunlight causes photoelectrons to be 
ejected. For most substances, the photoelectron current is 
on the order of one to four nanoamps per square centimeter. 
Since this is comparable to or larger than most space plasma 
electron currents, the surface will tend to acquire a small 
positive charge. The actual equilibrium potential will 
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depend on the details of the ion and electron distribution 
function, however (Whipple, 1976). The fluxes to a sunlit 
plate immersed in a plasma are shown schematically in Fig. 

2. The lit side will tend to charge slightly positive, and 
the dark side negative. 

The NASA MSFC baseline design (Hanley, 1978) is shown 
in Fig. 3. The surfaces on the satellite are divided into 
two types: active and passive, depending on whether or not 

voltages appear on the surface as a result of the satellite's 
own power supply. Passive surfaces include the solar re- 
flectors and structural members. Active surfaces include 
the solar cells, interconnects, and bus bars. Active sur- 
faces may attract or repel the ambient ions or electrons 
depending on the polarity of the surface voltage. Currents 
reach the passive surfaces only by photoemission and the 
thermal motion of ions and electrons. (We ignore backscat- 
tered and secondary electrons.) 

Calculation of Potentials 


We make the simplifying assumption of a thin sheath (or 
1-dimensional) approximation, i.e., the area collecting 
plasma is the actual geometrical area of the satellite (no 
focussing considered). The ambient electron and ion cur- 
rents, therefore are, simply the thermal currents, given by 


'i,e 


N 


( 


8kT J/2 


ttM 


) 


( 1 ) 


where N, e, T, and M are the number density, charge, temper 
ature and mass for electrons or protons, depending on which 
current is calculated. 


Parker (1979) has addressed the problem of a large flat- 
plate solar collector in space. He has found that the thin- 
sheath approximation is not valid at geosynchronous orbit 
for active structures. However, in the MSFC design, the 
passive, grounded reflecting panels form a trough in which 
the solar cells lie. Since the reflectors are conducting, 
tne.y have a tendency to confine electric fields from the 
solar cells within the trough. This reduces the thick-sheath 
focussing effect because the electric fields do not pene- 
trate signif iciantly into space above the trough, and the 
reflectors themselves are barriers against plasma fluxes 
entering from the sides of the trough. Later in the paper 
we verify this assumption by showing results from a modified 
version of Parker's PANEL program for the special geometry 
of the MSFC design. 



The analytic calculations below assume, for simplicity. 
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an intermediate sheath approximation; I.e., no focussing of 
outside plasma is considered* yet the sheath is large enough 
that photoeleetrons from the reflectors can impact the solar 
cell, and vice versa. 

For GEO, our assumed "worst case" plasma conditions 
are: Nq = Ni ■ 2/cm 3 , kTe » 5 key and kTf ■ 10 keV (Inouve, 
19761. This yields Jq ■ 3 x 10‘ 10 A/cm 2 and Ji * 1 x 10” 11 
A/cmS 


The photoelectron current density was calculated by 
integrating the product of the photoelectron yield function 
for synthetic sapphire and the solar spectrum: the resul- 

ting photocurrent density Jpq is 3 x 10“ 9 A/cm. A similiar 
calculation for aluminum yields roughly the same photoelec- 
tron current density. 

It is apparent, then, that the photoelectron current 
will usually dominate for all sunlit surfaces at GEO. The 
•equiKbr+uni potential for such surfaces will be on the order 
of a few times the average photoelectron energy, from about 
1 to 100 V positive, such as is found on the dayside of the 
moon (Reasoner and Burke, 1972; Freeman and Ibrahim, 1975). 
Passive sunlit surfaces will attain this voltage; however, 
for active surfaces, the finite conductivity of the cover 
surfaces (kapton and sapphire) will prevent this voltage 
from being obtained, i.e., the surface potential will more 
nearly follow that of the underlying solar cell. 


Nightside potentials are estimated from Chopra's (1961) 
equation: 


.♦i - 7^ 1" ( 


HjTe ) 

M e Ti 


( 2 ) 


For the "worst case" described above, this implies a dark- 
side potential of -17,000 V. Secondary electron emission or 
backscattering will reduce this potential somewhat. Again, 
passive surfaces will attain this voltage, but most active 
surfaces will be more nearly the potential of the underlying 
solar cell . 


The most vulnerable active surfaces on the satellite 
are the solar cells because the ohmic contacts are separated 
from the plasma by only tens of micrometers of shielding. 
Figure 4 shows the dimensions and structure of the solar 
cell selected in the MSFC design. The GaAlAs cell is sup- 
ported from below by a kapton blanket and is covered with 
synthetic sapphire. The sapphire coverglass is 20 urn thick 
and the kapton blanket is 25 vm thick. 



For our study, the solar cell was idealized as a sap- 
phire - active region - kapton sandwich as shown in Fig. 5. 
Plasma ions were assumed to be attracted to the negatively 
biased porton of the solar array and plasma electron.s to the 
positively biased portion. Photoelectrons were assumed to 
leave the negative surface and be attracted to the positive 
surface. Secondaries were neglected. The currents used 
were those described previously; we assume a steady state 
condition. In this case the voltages across the sapphire 
and kapton dielectrics are the photoelectron and plasma cur- 
rents multiplied by the resistance pf the dielectrics. The 
assumed resistivity of sapphire is 10 1Z ohm-cm. Based on 
the measurements of Kennerud (19'/4) we have approximated the 
resistivity of kapton by 

p 8 9.2 x 10 16 exp -[E/1.1 KV/mil] ohm-cm, 

where E is the electric field across the kapton in KV/mil. 
The transcendental equation for the potential difference, V, 
through the 1-mil kapton layer is In [V/K] = -V/1100, where 
K is proportional to tne current (K = 9 xlC ife x thickness 
(cm) x current (A/cm ). This equation was solved numeri- 
cally. The resulting voltages are shown on Fig. 5: a drop of 
949 V through the ion-attracting side, and a drop of 3.3 KV 
through the electron-attracting side. In no case are the 
breakdown voltages exceeded; however, the voltage on the 
positive 8 array is within a factor of 2 of the breakdown volt- 
age. For an electron current ten times larger (which can 
certainly occur within the satellite's life-span), the volt- 
age drop is 5.4 kV, which is near breakdown. For this 
reason, we recommend replacing kapton with a higher conduc- 
tivity material, or else providing a current path from the 
solar cell to the back side. Conductive coatings will also 
help reduce spot arcing (McCoy and Konradi, 1979). 

Kennerud (1974) and others have found anomalous arcing 
when solar panels are held at high voltage negative in a 
plasma. Typical voltages and currents required for such 
anomalous arcing to take place are 400 volts at 1 x 10" 7 
A/cm 2 . Our expected ion currents_^o the negative portion of 
the solar array at GEO are 1 x 10 7 1 A/ cm 2 . Therefore, we 
do not anticipate anomalous arcing in the GEO environment. 

The MSFC design calls for the reflectors to be con- 0 
structed from 0.5 mil (12.5 urn) kapton covered with a 400 A 
film of aluminum. We expect the aluminized front side 
potential to be fixed at 1 to 100 volts positive by photo- 
electron emission. Using the same analysis that was applied 
to the kapton solar cell blanket, we calculate the reflector 
back side voltage to be approximately -1.7 kV for our stan- 
dard "worst case" condition, and-2.7 kV for a ten times 
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larger electron current. The breakdown voltage for 1/2 mil kapton Is 3.1 kV, 
which could be reached with only slightly more severe plasma conditions. 

Cl.arly, the backside must also be conducting and electrically connected to 
the front, or the kapton must be replaced with a higher conductivity material. 

A summary of the expected voltages on various surfaces during sunlit and eclipse 
conditions is shown on Fig. 6. Note that during eclipse the entire satellite 
may charge to high voltage negative. This should be countered by the use of a 
hot filament electron emitter to bleed electrons from the spacecraft. 

Optimizing the Grounding Point 

The currents between the satellite and the plasma will 
adjust until the net current is zero. This means that the 
flow of current to the positively biased areas must equal 
that from the negatively biased areas. In a flate plate 
collector, the balance is between plasma electron currents 
to the positive portions and plasma ion currents to the 
negative portions of the array. Since the plasma electron 
currents are so large, the plate will "float" substatially 
negative, thel-area of the collector with negative 

potential is much larger than the corresponding positive 
potential area (Parker, 1979). 

In the MSFC design, however, the large aluminum reflec- 
tors are also sources and sinks of photoelectrons. Photo- 
electrons from the reflectors will be attacted to positive 
portions of the solar cell array and photoelectrons from the 
negative portions of the solar cell array will be attracted 
to the neighboring reflector (Figure 7). These electrons 
will "hop"’ along the surface (Pelizarri and Criswell, 1978), 
adding to the power drain. Thus the photoelectron current 
becomes the dominant parasitic current, at least in all but 
the most intense substorm environments. 

The large aluminum reflectors make a convenient space- 
craft ground, since the sunlit sides will remain a few volts 
positive with respect to space. To minimize the power drain, 
the solar cell array should drive no new currents through 
the reflectors to the plasma. Thus the reflector "ground" 
should be tied to the solar cells in an optimium way. Accu- 
rate calculation of the 3-dimensional electric field pattern 
and resultant power drain including effects of the space 
charge and secondaries is a formidable task; an oversimplifed 
argument follows. If A- is the solar cell area that is 
negative and A+ is the solar cell area that is positive, 
current balance requires 

(A-) (J pe ♦ 20, ) - (ft*) <J pe + *J e > < 3 > 

A-/A+ = (Ope + 2J e )/(Jp e + 2J^ ) . 


or. 



Here we assume that the photoelectron flux from the reflec- 
tors to the positive segments Is approximately the same as 
the photoelectron flux from the negative segments to the 
reflectors. For low plasma-current regions, (e.g., the 
plasmasphere or the quiet plasmasheet) or for cases in which 
the plasma current Is shleloed from the surfaces magneti- 
cally, the ratio approaches unity. Even for our "worst 
case," the ratio Is only 1.17. Therefore, we recommend 
grounding the midpoint of the string to the reflectors. On 
the other hand, at low Earth orbit plasma electron and ion 
ram fluxes dominate, and the grounding point must be more 
carefully calculated. 

With the ground point determined, the parasitic load 
can be calculated. The principal parasitic current at GEO 
is from photoelectrons (Fig. 7), and is calculated to be 
about 3000 A . Coupled with an average potential drop of 
11375 V, this implies a power loss of 34 MW, which is only 
0.7% of output power, and is easily manageable by slight 
oversizim,. This percentage power loss is comparable to 
that ( * 0.1%) from a flat-plate collector (Parker, 1979). 
Thus optimizing the grounding point at GEO is not critical. 
As discussed later, however, at LEO optimization could be 
very important. 

Currents at Low-Earth Orbit 


An integral part of the SPS concept is the Earth-Orbit 
Transfer Vehicle (E0TV) which will transfer the SPS to GEO. 

It is expected to employ a high-voltage solar cell array and 
to operate primarily in the low-Earth orbit (LEO) environ- 
ment where the plasma currents are considerably different 
than GEO. At 400 km altitude, the dominant ion is 0 + with a 
number density of 10 6 /cm 3 and a temperature of 2000 K 
(Johnson, 1965). Thus the thermal ion current will be 7 x 
10" 9 A/ cm and the thermal electron current will be 3 x 10" 7 
A/cm . For these currents, the potentials on the EOTV will 
be comparable those for which Kennerud (1974) found arcingjthere- 
fore, one must expect arcing to take place on negatively- 
biased surfaces unless a lower-voltage array is used. 

Indeed, arcing has been observed from insulated surfaces in 
a LEO simulation vacuum tank test (McCoy and Konradi, 1979). 
Alternatively, the array could be biased with a minimum of 
negative surface (grounding the lowest end of the string to 
the reflectors), but that would be far from the optimum 
grounding scheme, and would increase parasitic losses by a 
factor of three. 

Spacecraft motion implies a substantial though varying 
ram flux which will cause an additional parasitic current 
drain of as much as 2 x 10“ 7 A/cm 2 . Coupled with the cur- 
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rent losses due to the thermal currents, the power loss could be as 
high as -3%. As noted, however, arcing probably will occur at much 
lower potentials than those for which 3% power loss would be observed. 

Parker (1979) has pointed out that sheath and wake effects also could 
substantially alter the satellite potentials and current flow. 

EOTV Parasitic Load Due to Thurster Charge-exchange Ions 

An additional source of parasitic cur-ent for the EOTV Is created by 
charge-exchange of Ionized neutral gas from the thrusters with the 
energetic Ions from the main thruster Ion beam. This results In "thermal" 

Ions which may drift Into the Langmuir sheath electric field region 
surrounding the solar cell array. Once Into the field they will be 
accelerated toward the solar cells and produce a parasitic load. 

Following an approach outlined by H. R. Kaufman (NASA CR-135099) we 
have estimated the resulting parasitic load to the EOTV solar array to 
be 174 MW or 52% (Freeman and Few, 1979). This load Is clearly Inacceptable 
but It can easily be mitigated by placing a shield between the thrusters 
and the solar array. This shield can consist of an alumnlzed kapton sheet 
stretched across the end of the EOTV. The shield will need to have a 
height comparable to the dimensions of the Langmuir sheath, about 500 m. 
Additionally the low voltage edge of the solar cell array should be located 
toward the outside. Similar shields should be considered adjacent to the 
ACS thrusters on the SPS itself. 

Non-Steady State 

Until now it has been assumed that the charging currents from the plasma 
are steady. This approach Is supported by a study of the time dependent 
charging of a three-axis stabilized spacecraft by Massaro et al . , (1977). 

For all the surfaces modeled, they found that the greatest differential 
voltages occurred In the steady state limit, although nearly Instantaneous 
changes in absolute potential were observed. However, In order to evaluate 
the effects of nori-steady charging, we calculated the RC time constant or 
discharge time of the relevant insulators, sapphire and kapton. The RC 
decay time is pc where p Is the resistivity and e the permittivity. For 
kapton this implies a time constant of 1 hr; for sapphire, 1 sec. Large 
magnetospheric changes can occur with 1 min - 1 hr time constants (Mcllwain, 

1974; Inouye, 1976). Therefore, high voltages can build up on the kapton in 
time intervals short compared to the discharge time. Transient charging is 
not expected to cause differential charging in excess of the steady state 
predictions, nevertheless, the large kapton time constant reinforces the pre- 
vious conclusion that kapton should be replaced with a higher conductivity material. 

3-Dimensional Model 

All of the foregoing analysis on parasitic loads, plasma induced voltages, 
etc., is based on one-dimensional plasma theory. More precise results require 
a three-dimensional self-consistent computer model which takes into account all 
plasma sources and interactions with reflectors simultaneously. A computer 
program, "PANEL"" written by Dr. Lee Parker (Parker, 1979), provided a convenient 
starting point for our model of the SPS environment. Preliminary results will 
be presented here. They are preliminary since we have not yet included the photo- 
electron current (which we showed to be important), nor have we as yet included 
space charge effects. Nevertheless, the results demonstrate several important 
features of the sheath around the SPS troughs. 



PANEL utilizes a three-dimensional grid wiierr the sat- 
ellite Is modeled by fixing potentials at selected grid 
points. Laplace's equation Is then satisfied by relaxing 
the free space potentials until Gauss's law Is satisfied for 
a box surrounding each point. The currents and power losses 
ire obtained by numerically performing the integral 


r 

0 


dv 


n/z 

/ do 
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2lf 


d$ f (v,6,4>) v 3 cose sine 


where J Is the current density, and f is the distribution 
function. The problem is then to evaluate f. For a colli- 
sionless steady state plasma, the Vlasov equation 

v * vf + ■“ f • $ y f ■ 0, states that a distribution function 

is constant along a particle's path in phase space. If f is 
written In terms of a particle's total energy (E ■ T + V, 
the kinetic plus potential energy), f will be constant in E 
along the path in real space* The Integral for J is then 
transformed into a sum using the method of gauss Ian gradu- 
atures which picks key values of E, 6, and 4>. These values 
represent trajectories that are traced backwards to either 
source or nonsource regions to determine the value of f. 

Once the current is known It is multiplied by the local 
potential to determine the power loss at that point. 

PANEL is a Laplacian calculation since space charge 
effects are not included in the electrostatic potential 
calculation. The next phase in the development of PANEL is 
to calculate the charge density for each point in space by 
evaluating the integral 


N * / dv / de / d 4»f ( v , 0,<}>)v 2 sine 
ooo 

in the same manner as described for the current calculation. 
Then PANEL must iterate between the potential relaxation 
routine and the density calculation since the density calcu- 
lation depends upon the potential structure for accurate 
trajectories. This is known as the inside-out method 
(Parker, 1977) because trajectories are traced backwards in 
time. 


Figure 8 illustrates the three dimensional grid used to 
model two interior panels of a trough. Not shown are grid 
points at the intersection of all integer x and y values and 
even values of z. One unit of grid sparing corresponds to 
85.0 meters, giving model dimensions of 765 m X 4I?5 m. Fixed 
voltages are indicated on t.be fiyure. The assumed plasma 



conditions are N< ■ N e ■ 2/cm 3 , kT ^ ■ 10 keV, kT e ■ 5 keV. 
For these conditions, the random thermal current densities 
are, as before: 

J th|1 - 1.25 x 10“’ A/m 2 

J t h,e ■ 3.79 * 10" 6 A/m 2 

The dimensionless numbers at selected points on the panels 
are ratios of local average electron current densities to 
the random electron thermal current. For the two panels 
modeled, PANEL traced 864 trajectories per grid square of 
surface. The resulting total current collected and power 
loss are 6.64 x 10'\A, and 5.66 x 10 2 W for protons and 
2.25 A and 2.72 x 10 4 W for electrons. Calculated potential 
patterns In the x ■ 0 plane and y ■ 3 plane are shown In 
Figs. 9 and 10, respectively. Note that potentials of only 
1-2 kV extend beyond the upper limits of the trough, justi- 
fying our earlier "intermediate sheath" approximation. 

Photoelectrons from the reflectors and backscattered 
and secondary electrons undoubtably will be important con- 
tributor to the power loss but have not yet been modeled. 

Magnetic protection of the SPS 

. The SPS of necessity contains bus bars of current 
KT A, routed between the solar panels and the microwave 
antennae. With judicious routing of these bus bars, the SPS 
can create its own protective magnetic barrier, screening 
out all the low energy (MOO eV) plasmaspher ic plasma (which 
can cause power drain), and most of the energetic electrons. 
Parker and Oran (1979) have shown that this idea is feasible 
with nominal bus-bar currents. We propose modified bus-bar 
currents to prevent spacecraft fields from merging with the 
earth's magnetic field. Merging can have two harmful ef- 
fects: 

‘0 It can channel energetic particles trapped in the 
Earth's magnetic field towards sensitive areas of the SPS. 

2) It can energize the high density plasmaspheric 
plasma that would otherwise be harmless. 

Previous spacecraft were small in size compared to 
particle gyroradii, so magnetic effects were not important. 
The size of the SPS, however, is comparable to particle 
gyroradii, so magnetic effects must, be taken into account. 
(At geosynchronous orbit, a 2 eV proton or 3 keV electron 
has a gyroradius of 2 km; a 50 eV proton or 80 keV electron 
has a gyroradius of 10 km.) In the following, in order to 
estimate these effects (i.e., to repel trapped particles and 
to minimize energy released in magnetic merging) we assume 



that it Is important to have spacecraft magnetic fields 
parallel to sensitive areas (e.g., solar cells) and aligned 
with the Eanh's magnetic field. (Even magnetic fields 
perpendicular to the surface can be beneficial, however, and 
have been considered In Parker and Oran, 1979). 

A solenoidal bus-bar winding yields the best magnetic 
field configuration: at a distance, the field approaches 

that of a dipole, and In the vicinity of the satellite the 
field is parallel to the solar panels. The windings tor the 
solenoid should enclose as much area as feasible. This will 
have two benefits: It will maximize the overall dipole 

moment while minimizing the bus bar length and thus 1R 
losses, and will minimize the internal field. On the other 
hand, for spatial uniformity, one should have a least one 
turn per kilometer. Some possible cross-sections are shown 
in Fig. 11. This figure Is a view from the north end of 
three types of trough- like SPS design and shows one turn of 
the helical winding each. 

The field of the SPS must have sufficient rigidity to 
successfully deflect the species desired to be excluded. 

Table 1 show magnetic moments v required for various tasks. 
Two possible orientations of the SPS's dipole moment are 
compared: parallel or antiparallel to the Earth's dipole 
moment. A parallel orientation, since it adds to the local 
magnetic field, is more efficient at shielding the SPS from 
particle bombardment; however, the opposite orientation is 
dynamically more stable, since the SPS's moment will tend to 
align with the Earth's magnetic field. In fact, the moment 
m a y be used to balance gravity-gradient torques if the dipole 
moment is large enough. -For a (uniform) body 22 km long and 
4 km wide of mass 5 x 10 kg, the moment of inertia about an 
axis perpendicular to the length of the satellite would be 
2 x 10 1 * kg-nr . The daily + 10 deg tilt of the geosynchro- 
nous magnetic field would ca"use a torque on the satellite of 
( ij x B ) ° 1.7 x 10 4 Nt-m, for a u = 10 12 A-m 2 (corre- 
sponding to 0.9 Nt of force on each end). Since the sat- 
ellite is so massive, this torque will result in a daily 
sinusoidal tilting motion of the satellite of amplitude * 

10" s degrees, completely negligable. A 10 deg tilt of the 
satellite toward the Earth, in contrast, will cause a grav- 
ity-gradient torque of 2.7 x 10 6 Nt-m, or 125 Nt at each 
end, requiring a magnetic moment of 1.5 x 10 14 A-m 2 to bal- 
ance it. Then, however, the 10 deg misalignment between the 
spin axis and the dipole axis of the Earth would become more 
important. In addition, the magnetic fields in the SPS 
center would be quite large (90 G.). The internal field is 
sensitive to the exact configuration, and can vary by a 
factor of two or sp depending on the area and number of turns 
per km. The rigidity, on the other hand, is not too sensi- 



tive on. the exact configuration, being mainly a function of 
overall magnetic moment. 

One reasonable magnetic field configuration Is shown In 
Figs. 12 *• 14. The dipole moment assumed for these figures 
Is the low-field case, 10 11 A-m 2 per km, 21 km total. All 
components of the field are, of course, linear In the dipole 
moment. This model superposes 21 dipoles at 1 km Intervals 
(simulating one turn per km). Figure 12 shows vector mag- 
netic fields for one quadrant; Fig. 13 shows contours of 
constant |B|» and Fig. 14 shows magnetic field components. 
Here the z-component is measured along the long axis and 
the p component is measured from the long axis. The center 
of the SPS is the lower left corner (z ■ 0, p * 0). Only 
one quadrant is shown because of symmetry: B z (z) * B z(-z); 
Bp(z) ■ -Bp(-z). The field is similar to that of a solenoid 
and is nearly parallel to the long sides of the SPS (and 
therefore to the solar cells), converging at the SPS's north 
and south ends. (The SPS is aligned north-south to minimize 
Jtbe shadowing of one SPS on another in the equinox seasons.) 

The field in Figs. 12 - 14 is strongest at the ends and 
weakest in the center; therefore, fewer wraps (or, more 
likely, less current per wrap) could be used at the ends and 
still obtain the same overall rigidity. A field of 100 
extends to over 7 km from the center, and a field of 20 
extends to 19 km. The overall rigidity at p = 1 km, z = 0 
km is roughly 2000Y - km (G-cm). With a magnetic field of 
this orientation and strength, ions < '10 eV (including all 
the plasmaspheric plasma) and electrons < 30 keV (most of 
the plasma sheet electron fluxes) are excluded. Higher 
dipole moments would yield more shielding (see Table 1). 
Thus, it appears that magnetic protection is feasible. 
Because of the convergence of the field, particle fluxes 
will have a tendency to strike only the ends of the long 
axis of the SPS. Simply capping the ends of the SPS, then, 
will be sufficient to protect electronics and humans inside 
from the lower-energy particles. Such capping is also useful 
to prevent the plasma from the ion engines from returning to 
the satellite, causing a significant power drain (Freeman 
and Few, 1979). 


Conclusions 

The SPS will certainly interact with its plasma envi- 
ronment. It appears that, with relatively minor modifica- 
tions to the NASA MSFC baseline design, these interactions 
will not significantly impair SPS operations. The conclu- 
sions and recommendations of this study include: 

1) Arcing is likely to occur on kapton surfaces (the 
solar reflectors and the solar cell back surface blanket) 
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during substorms unless the kapton Is replaced by a lower resistivity 
material (p <_ 10 13 ohm-cm) or current paths from the surfaces to the 
solar cells are provided. 

2) The SPS parasitic load under normal conditions will be about 34 MW 
(for a 5 GW array) at geosynchronous orbit. This 0.7% power loss should 
be accomodated by oversizing. 

3) The optimum grounding point at GEO for the SPS solar cell array Is 
approximately the midpoint on the voltage string. At LEO, arcing considera- 
tions demand that the string be biased mostly positive, although the optimum 
configuration to minimize power loss would be substantially negative (see 
conclusion 5). 

4) The solar cells may require conductive coatings. The reflector 
panels may require current paths linking the front and back sides. Laboratory 
tests in a substorm simulator on realistic solar panels are recommended to 
determine the actual arcing probability. 

5) Severe arcing problems are expected for negative portions of the EOTV 
solar cell array at LEO. Overcoming this problem by biasing the array as 
positive as possible will result in high parasitic loads (power losses on 
the order of 3%). Only a low voltage EOTV solar array should be used. 

6) The SPS will occasionally charge to about -20 kV during eclipses. An 
active discharge method such as a hot filament electron emitter should be 
provided. 

7) A shield should be placed across the ends of the EOTV to prevent 
thruster ion feedback to the solar array. Similar shields may be required 
on the SPS. 

8) Three-dimensional computer modeling of the SPS electric field pattern 
and plasma currents is underway. The model shows that, for the grounding 
scheme used here, spacecraft electric fields extend only slightly beyond the 
reflectors. 

9) Active magnetic plasma shielding is possible through judicious routing 
of bus-bars; power drain from additional lengths of bus-bars has not been 
calculated yet. 

10) It is possible to use the internal magnetic field to align the satellite 
(counteracting gravity-gradient torques), but it would require an unreasonably 
large magnetic moment (1.5 x 10 14 A-m 2 ). 

Acknowledgments 

The authors thank Dr. Lee Parker for consultation and the use of the 
computer program "PANEL." In addition, we have benefited from discussion 
with Dr. James McCoy. This work was supported by NASA under grant NAS8- 
33023. 

The computer program PANEL is attached as appendix A. 



References 


Burke, W. J., Reiff, P. H., and Reasoner, D. L., "The 
Effect of Local Magnetic Fields on the Lunar Photoelectron 
Layer While the Moon Is in the Plasma Sheet," Geochemicha 
Cosmochemicha Acta, Suppl. 6 , Vol. 3, 1975, pp. 2985-2997. 

Chopra, K. P., Reviews of Modern Physics, Vol. 33, 1961, p. 
153. 

DeF orest, S. E., "Spacecraft Charging at Synchronous 
Orbit," Journal of Geophysical Research, Vol. 77, February, 
1972, pp. 651-659. 

Freeman, J. W. and Few, A. A., "Electrostatic Protection of 
the Solar Power Satellite And Rectenna," final report, NASA 
contract NAS8-33023, Rice University, Houston, Texas, May 1979. 

Freeman, J. W. and Ibrahim, M., "Lunar Electric Fields, 
Surface Potential And Associated Plasma Sheets, The Moon, 

Vol. 14, 1975, pp. 103-114. 

Garrett, H. B. and DeForest S. E., "Analytic Simulation of 
the Geosynchronous Plasma Environment," Planetary Space 
Science, Vol 26, (in press), 1979. 

Hanley, G. M., "Satellite Power Systems (SPS) Concept Def- 
inition Study," final report for NASA contract NAS8-32475, 
Rockwell International Report #5078-AP-0023, Cincinnati, 

1978. 

Inouye, G. T., "Spacecraft Potentials in A Substorm Envi- 
ronment," Progress in Astronautics and Aeronautics ; Space - 
craft Charging by~~Magnetospher ic Plasmas , ed. by'A. Rosen, 

A1AA, New York, 1976^ pp, 103-120. 

Johnson, F. S (ed). Satellite Environment Handbook , Stan- 
ford University Press', Stanford, California, 1965. " 

Kennerud, K. L., "Final Report High Voltage Solar Array 
Experiments," Boeing Aerospace Corp., Seattle, Wash., Rept. 

No. CR121280, 1974. 

Massaro, M. J., Green, T., and Ling, P., "A Charging Model 
for Three-axis Stabilized Spacecraft," Proceedings of the 
Spacecraft Charqinq Technology Conference, edited by C. P. 

■pTke and R. R. Lovell, Air Force Rept. AfGL-TR-77-0651 , 1977, 
pp. 237-267. 



McCoy, J. E. and Konradl, A., "Sheath Effects Observed on 
A 10 Meter High Voltage Panel In Simulated Low Earth Orbit 
Plasma," Space Craft Charqinq Technoloqy-1978, edited by R. 

C. Finke 'and'T.' P Pi kcC NAOTP 2071, l$737“pp. 315-340. 

Mcllwain, C. E., "Substorm Injction Boundaries," Magneto - 
spheric Physics , edited by B. M. McCo.. r.-c, D. Re i del, 

ITingham, Mass., 1974, p. 143. 

Mizera, P. F. and Fennell, 0. F., "Satellite Observations 
of Polar, Magnetotail Lobe, and Interplanetary Electrons at 
Low Energies," SAMSO Rept TR-78-6, 1978. 

Parker, L. W., "Plasma Sheath Effects and Voltage Distri- 
butions of Large High-power Satellite Solar Arrays," Space - 
craft Charginq Technology-1978, edited by R. C. Finke and C. 
P. Pike, NASA CP 2071, 1979, pp . 341-357. 

Parker, L. W. and Oran, W. A., "Magnetic Shielding of Large 
High-power-satellite Solar Arrays Using Internal Currents," 
Spacecraft Charqing Technology-1978, edited by R. C. Finke 
and C. P. Pike, NASA tP 2071, 1979, pp. 376-387. 

Parker, L. W., "Calculation of Sheath and Wake Structure 
about A Pillbox-shaped Spacecraft in A Flowing Plasma," 
Proceedings of the Spacecraft Charging Technology Conference , 
edited by C. P. PTke and R. R. Lovtll, Air Force Report 
AF6L-TR-77-0051 , 1977, pp. 331-336. 

Pelizzari, M. A. and Criswell, D. R., "Differential Photo- 
electric Charging of Nonconducting "Surfaces in Space," Jour- 
nal of Geophysical Research, Vol. 83, November, 1978, pp. 
5233-5244. 

Pike, C. P. and Brun, M. H., "A Correlation Study Relating 
Spacecraft Anomalies to Environmental Data," Progress in 
Astronautics and Aeronautics : Spacecraft Charging by Mag - 

netospheric Plasmas, edited by A. Rosen, AIAA (New York), 
T976, pp. 45-60. 

Reasoner, D. L. and Burke, W. J., "Characteristics of the 
Lunar Photoelectron Layer," Journal of Geophysical Research , 
Vol. 77, December, 1972, pp. 6671. 

Reasoner, D. L., Lennartsson, W., and Chappell, C. R., 
"Relationship Between ATS-6 Spacecraft-charging Occurences 
and Warm Plasma Encounters," Progress in Astronautics and 
Aeronautics ; Sapcecraft Charging by Magnetospheric Plasmas , 
edited by A. Rosen, AIAA, New York, 1976, pp. 89-101. 



Reiff, P. H., "Magnetic Shadowing of Charged Particles by 
An Extended Surface," Journal of Geophysical Research, Vol. 
81, July, 1976, pp. 3423-3427; 

Reiff, P. H. and Burke, W. J., "Interactions of the Plasma 
Sheet with the Lunar Surface at the Apollo 14 Site," Journal 
of Geophysical Research , Vol 81, September, 1976, pp. 4761- 
4764. 

Shaw, R. R., Nanevicz, J. W., and Adamo, R. C., "Observa- 
tions of Electrical Discharges Caused by Differential Charg- 
ing," Progress in Astronautics and Aeronautics : Spacecraft 

Charging by Magnetospheric Plasmas, edited by A. 'Rosen, AIAA, 
New York, 1976, pp. 61-79T 

Whipple, E. C., Jr., "The Equilibrium Electric Potential 
at A Body in the Upper Atmosphere and in Interplanetary 
Space," NASA T. N. X-615-65-296, 1965. 

Whipple, E. C., Jr., "Theory of the Spherically Symmetric • 
'Photoelectron Sheath: A Thick Sheath Approximation and 

Comparison with ATS6 Observation of a Potential Barrier," 
Journal of Geophysical Research, Vol. 81, February, 1976, 
pp. 601-607. 



FIGURE CAPTIONS 

Fig. 1 Sketch of the Earth's magnetosphere (from Mizera and 
Fennell, 1978). 

Fig. 2 Schematic of plasma 
and photoelectron currents. 

Fig. 3 Sketch of the MSFC January 25, 1978, baseline design 
(from Hanley, 1978). 

Fig. 4 Cross-section of a proposed GaAlAs solar cell (from 
Hanley, 1978). 

Fig. 5 Idealization of the solar cell blanket, used in cal- 
culations of electrostatic potential, for the "worst case" 
plasma fluxes. 

Fig. 6 Summary of voltages on the reflectors and solar cells 
surfaces, for solar cells at large positive voltages (top), 
large negative voltages (miodle), and ouring eclipse (bot- 
tom). (Midpoint of the solar cell voltage string is assumed 
to be grounded to the sunlit side of the reflectors.) 

Fig. 7 Summary of parasitic current densities for the SPS 
and the parasitic current and power loss total for one half 
of the Marshall satellite (5 GW system). 

Fig. 8 Computer grid used to model 2 panels of the SPS. 
(Small numbers on the panel surface are the plasma electron 
currents normalized to random thermal currents.) 

Fig. 9 Equipotential contours in the yz plane at x * 0 (in- 
dicated in Fig. 8). 

Fig. 10 Equipotential contours in the xz plane at y = 3 (in- 
dicated in Fig. 8). 

Fig. 11 Recommended current windings for several SPS con- 
figurations (view from north end). 

Fig. 12 Vector magnetic fields 
for a solenoidal current 
configuration, low-field case 
( u = 10 11 A-m 2 per km, 21 km 
total). (Z-axis is along the 
spacecraft (z = 0 is the 
center), and p is measured 
from the spacecraft axis; only 
one quadrant is shown, because 
of symmetry: Bz(-z) = Bz(z); 

B p (-z) = -B p (z) . ) 
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Fig. 13 Contours of constant 

i Bj for the low-field case, 
niy one quadrant Is shown. 

Fig. 14 Contours of constant B p and B z » low-field case. 
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Magnetic Moment Required for SPS Tasks 


Task 

Rigidity 

Required 

Orientation 
Of Moment 

Internal 

Field 

(Gauss) 

Required 

Moment 1 

A-mVkmt 

Shielding 
200 eV 
Protons 
and 30 keV 
Electrons 

2 x 10 3 

^Parallel 

Antiparallel 

1.3 

4 

1 x 10}} 
3 x 10 n 

Shielding 
3 KeV 
Protons 
and 2 MeV 
Electrons 

8 x 10 3 

^Parallel 

*Antiparallel 

5.3 

11 

4 x 10 11 
8 x 10 11 

.Shielding 
30 KeV 
Protons 
10 MeV 
Electrons 

3 x 10 4 

Parallel 

*Antiparallel 

20 

25 

1.5 x 10 12 
2 x 10 12 

Magnetic 
Alignment 
(Bai ance 
Gravity- 
Gradient) 

N/A 

*Anti parallel 

92 

7 x TO 12 


Recommended Orientation 

■fMultiply by 21 for total magnetic moment. 
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VOLTAGES SHOWN ARE RELATIVE TO THE LOCAL SOLAR 
CELL VOLTAGE. THE',’ REPRESENT THE IR DROP ACROSS THE. 
COVER GLASS OR KAPTON BLANKET. 

• THE KAPTON BREAKDOWN VOLTAGE IS * 6250 V 
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Figure 6 




IH£-VARI_QU$ XURRENT DENSITIES ARE : 


J PHE = 3x 1CT 9 amp/cm 2 (FOR SAPPHIRE) 
J e CPLASf^A) = 3 x 10’ 10 amp/ cm 2 
(PLASMA) -lx 10”^" amp/cm 2 



TOTAL PARASITIC CURRENT : 
*p = 3000 AMPS 
V = 11,375 V 


FTOTDELFCmON FLO^f DIRECTION 


THE PARASITIC POWER IS : 

P p = 34 MW 
(0.7% OF OUTPUT POWER) 
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EQU I POTENTIALS IN THE Y = 3 PLANE 
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EQU I POTENTIALS IN THE X = 0 PLANE 
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PREFERRED CONFIGURATION 




AREA= 1.6 X 10 6 M 2 

LOOP LENGTH = 5200 M 

1-4 TURNS PER KM RECOMMENDED 



AREA = 9.1 X 10 s M 2 AREA= 2.4 X 10 G M 2 

LOOP LENGTH = 4300 M LOOP LENGTH = 9500 M 

1-6 TURNS PER KM RECOMMENDED 1-2.5 TURNS PER KM RECOMMENDED 


Figure 11 
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SUMMARY AND CONCLUSIONS 


1. The very high lightning flash density In :;.any parts of the United 
States and the large size of the SPS rectenna require us to Incorporate 
lightning protection systems In the rectenna design. 

2. A distributed lightning protection system Is described In this report 
that will protect the rectenna components from direct lightning strike 
damage and will. In addition, provide reduced Induced lightning effects 
In the power and control circuits. 

3. The proposed lightning protection system should be Incorporated as a 
structural member of the rectenna support system; viewed as such, the 
lightning protection system will not appreciably Increase the total 
material requirements for the rectenna unless materials are used that are 
Incapable of safely conducting lightning currents. 

4. The lightning protection design places the conducting elements so that 
the microwave shadow cast by protection systems falls along the upper 
edge of the billboard on which it Is mounted (and the lower edge of the 
next billboard to the north); these shadow areas are only a slight 
fraction of the collecting area, so the protection elements produce very 
little, If any, additional power loss to the rectenna as a whole. 

5. Individually the microwave diodes are self- protecting with respect to 
"average" lightning and those near the center of the rectenna are safe 
from extreme lightning. However, the series connection of the diodes to 
form 40,000 V strings creates a protection requirement for the string. 
Standard surge protection practices are necessary for the string. 

6. Electric power industries usually attribute 10% of the cost of power 
transmission equipment to lightning protection requirements. If this 
factor Is not already included In cost estimates, It should be added. 
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Part II 


LIGHTNING PROTECTION OF THE RECTENNA 
ABSTRACT 


Computer simulations and laboratory tests were used to evaluate the 
hazard posed by lightning flashes to ground on the SPS rectenna and to make 
recommendations on a lightning protection system for the rectenna. The 
distribution of lightning over the lower 48 of the continental United States 
was determined, as were the Interactions of lightning with the rectenna and 
the modes In which those Interactions could damage the rectenna. The studies 
showed that lightning protection was both required and feasible. Several 
systems of lightning protection were considered and evaluated. These included 
two systems that employed lightning rods of different lengths and placed on 
top of the rectenna s billboards and a third, distributed system. The distri- 
buted system Is similar to one used by power distribution companies; It con- 
sists of short lightning rods all along the length of each billboard that are 
connected by a horizontal wire above the billboard. The system that not only 
affords greater protection than the others considered but offers easiest 
Integration Into the rectenna' s structural design, was the distributed light- 
ning protection system. 
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SUMMARY OF THE RECOMMENDED LIGHTNING PROTECTION DESIGN 


Based upon our research, computer simulations, and laboratory tests with 
a scale model, we recommend a distributed lightning protection system that 
employs a horizontal conducting member with points and grounds placed at every 
bay or billboard (14.69 meters apart). This configuration not only provides 
greater protection than other configurations that were evaluated, It Is more 
easily Integrated Into the structural design of the rectenna. The recommended 
system Is shown In Figure 1. 
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DISTRIBUTED LIGHTNING PROTECTION SYSTEM 

FIGURE 
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PREFACE 


The objectives of this study are to evaluate the hazard posed by light- 
nlng flashes to ground on the SPS rectenna and to make recommendations for a 
lightning protection system that will provide sufficient protection to the 
rectenna. For purposes of this study, the SPS rectenna design Is based upon 
the data supplied to us by Rockwell International In July. 1978. 

This study has four major components, each with several elements of 
Investigation. The components were: lightning distribution; lightning 
Interactions; rectenna damage estimates; rectenna protection. The elements of 
each component are listed In Table A. The study plan was to proceed from top 
to bottom evaluating the elements listed In each component; work proceeded In 
a parallel manner for the four components. The organization of this final 
report reverses this order by presenting the more Important results of the 
study first, then following this with the material and considerations leading 
to the conclusions. 



TABLE A 

Rectenna Electrostatic Protection 


Lightning 

Distribution 

1. Obtain climat- 
ological data. 


2. Format data 
for computer 
use. 


3. Construct 
program for 
computation 
of lightning 
density. 

4. Produce contour 
map of light- 
ning density. 


Lightning 

Rectenna Damage 

Rectenna 

Interactions 


Estimates 

Protection 

1. Review/compile 

1. 

Diode failure 

1. Panel 

data on 


modes (scaled 

transient 

lightning 

parameters. 


from avail- 
able diodes.) 

protectors 

2. Construct 

2. 

Insulation 

2. Billboard 

program for 
computation 
of fields and 
currents in 
the rectenna 
plane from 
parameterized 
lightning 


breakdown. 

surge 

protectors 

3. Evaluate en- 

3. 

Direct strike 

3. Lightning 

hancement 

factors. 


damage estimates. 

Inverter 

protectors 


4. Conduct 

4. Direct strike 

4. Lightning 

1 aboratory 

damage estimates. 

rod 

simulations. 


systems. 
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I. ANALYSIS OF LIGHTNING ROD PROTECTION CAPABILITIES FOR A CONFIGURATION 
SUITABLE TO THE RECTENNA 

1. Cone of Protection Considerations: 

I. 1.1 Definition and Considerations 

The capability of a vertical conductor to attract a lightning flash 
Is described by the cone- of- protect Io n, or perhaps more accurately the 
cone-of-attractlon. In theory, any lightning flash that would have 
entered this cone had the vertical conductor not been In place, will 
strike Instead the conductor and be shunted to the ground. The method by 
which this process takes place Is as follows: 


The lightning stepped leader creates high voltages over a wide 
area on the rectenna because of the large charge on the leader 
tip. At points on the rectenna where the electric field reaches 
breakdown values due to local enhancement factors, upward 
propagating sparks are initiated which move to meet the downward 
propagating stepped leader. The upward propagating spark which 
first makes contact with the leader completes the electrical circuit 
and the lightning flash current will pass through the structure that 
initiated the successful upward going spark. 


The cone of protection Is primarily a function of the height of 
the vertical conductor because of the field-enhancement factor which 
enables the taller object to initiate the upward spark before lower 
objects. Other factors enter Into the consideration of the cone of 
protection, such as the charge on the leader tip and the velocity of the 
leader, because these factors strongly Influence the timing of the 
production of upward sparks and the height at which the spark and leader 
meet. In general, the results of research into this subject have shown 
that the larger the leader charge, then the larger the angle 8 of the 
associated cone of protection. Since larger leader charges are usually 
associated with the larger lightning currents, we find a fortunate result 
that the cone of protection increases with the potential hazard of the 
lightning flash. 

It follows then that the angle 8 of the cone of protection (See 
Figure 2) varies with the particular lightning flash. 8= 45° is a very 
commonly used design angle in the United States and many of the examples 
in this report employ 8 s 45°. 
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Figure 2 


1.2 Distributed Lightning Protection Systems 

The cone of protection and the experimental data used to evaluate 
are specifically related to the single elevated point, and in most cases 
the system under consideration is 10 to 100 meters in height. As will be 
seen later, lightning protection of the rectenna falls into a class of 
structures that requires distributed lightning protection tactics. 

Figure 3 illustrates a distributed system used by power transmission 
companies. The main point is that the cone of protection concept is of 
limited usefulness in the total protection problem. We will use it on 
the panel and billboard scale as a technique to make a comparative 
assessment of capabilities of various configurations. 

2. Lightning Rod Protection Configurations Compatible with the SPS Rectenna 

We have considered three different configurations of lightning rod 
systems in this effort. In the smallest scale system considered each rectenna 
panel (0.74m in width) had a short lightning rod attached; see upper example 
in Figure 4. In the medium scale system each rectenna support structure 
(14.69m apart) or billboard will have an attached lightning rod; see middle 
example in Figure 4. And, in the distributed protection system, short termi- 
nals located on each rectenna support structure (14.69m apart) were connected 
by horizontal conducting structures; see lower example in Figure 4. 

As seen in the analysis of the billboard scale system, it is impractical 
to seriously consider larger scale systems. 
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POWER LINES EMPLOY DISTRIBUTED LIGHTNING PROTECTION SYSTEMS. THIS ILLUSTRATION 
SHOWS A "STATIC" OR GROUNDED PROTECTION WIRE TAKING A STRIKE AND PROTECTING THE 
POWER LINES BELOW. 


FIGURE 3 
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2.1 Lightning Rod Protection at the Panel Seale 

In this system configuration a relatively short lightning rod Is 
positioned at the top of each panel and oriented perpendicular to the panel 
face (see Figure §). Conceptually the rod Is centered on the top of the 
panel, but In practice It could be on the panel edge without altering the 
results of this analysis. 

Let a be the Inclination of the rectenna. Figure 6 Illustrates the case 
where 8, the angle of the cone of protection, Is greater than o. This figure 
applies only to the conditions In the vertical plane that passes through the 
lightning rod and Is perpendicular to the rectenna face. In this particular 
projection it appears that the short (example 0.74m) lightning rod on the 
panel provides adequate protection to the rectenna. In other projections we 
see that there are, however, "holes In the armor. 

Figure 7 Is an enlargement (xlO) of the lightning rod portion of Figure 
6, and defines the parameters to be used In the following discussions. The 
cone of protection intersects the plane of the rectenna to form conic 
sections: 

(1) If a+ 90° the intersection is a parabola. 

(2) If o+ 0 < 90°the intersection is an ellipse. 

(this is the case Illustrated In Figures 6 & 7) 

(3) If a+ $> 90° the intersection Is a hyperbola. 


If we now look at the Intersection of the cone of protection with 
the panel for the specific cases above, we see the emergence of the 
protection problem with this type of lightning rod protection configura- 
tion. From the geometry of Figure 7 we see that the axis of the cone is 
at s. = L tan a and that the vertex of the conic Is at d s L 
tan ( 8 - a) relative to the top of the panel. 



Figure 6 



FIGURE 7 

ENLARGED VIEW OF THE UPPER END OF THE 
RECTENNA IN FIGURE 7 . 
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In a coordinate system defined In the reetenna plane with the origin 
at the axis of the cone and the y axis directed north (toward top of 
reetenna) and the x axis directed east, the equation for conic Is: 

2 2 2 “ " 


x cos a 
L 2 tan 2 



2 + 2 y sin a cos a D j 


L tan & 8 

For the parabolic solution this equation reduces to: 


2L sin 2 8 
cos B cos a 


(y- 


2 cos 8 COS a 


) 


In figure 8 we have plotted the Intersection of cones of protection for 
three lightning rods of lengths 0.185m (o V 4 panel width), 0.37m (° V 2 panel 
width), and 0.74m (° panel width.) 


In these examples the reetenna Inclination angle a Is taken to be 45° and 
the cone of protection B Is equal to 45°. The resulting Intersections are 
parabolas for the cases depicted In Figure 8 . For the parabolic solution the 
cone of protection Is parallel to the face of the reetenna In the vertical 
plane bisecting the panel (The view of Figure 6 and 7 except that 
here a ° B 0 45 u ) . 

At lower latitude sites (below 40°) the reetenna Inclination angle a is 
less than 45° and the 45° cone of protection intersection becomes an ellipse; 
In Figure 6 the vertical projection Illustrates the Intersection In the plane 
through the lightning rod. The elliptic solutions leave regions along the 
base of the reetenna unprotected. Henee, the parabolie solutions of Figure 8 
and the table (Fig. 9) represent maximum protection capabilities of the cone 
of protection with the panel scale protection configuration. The small 
ellipse In Figure 11 shows the cone of protection intersection 
for a = 40 B - 45 °, and L c 0. 74m. 

2.2 Lightning Rod Protection at the Bay or Billboard Scale 

In this system a longer lightning rod is placed at the center (or 
end) of each bay or billboard making them 14.69m apart. The mathematical 
description here is identical to that for the panel scale system ( 2 . 1 ). 
Only sizes are different. Figure 10 illustrates the billboard scale 
system. 


8 




FIGURE 8 

THE INTERSECTION OF THE CONE OF PROTECTION WITH A RECTENNA 
PANEL (THE CURVED LINES) SHOWN IN THE PLANE OF THE PANEL. 
LIGHTNING ROD LENGTHS = M, k AND 1 TIMES THE PANEL WIDTH 
ARE SHOWN PROJECTED VERTICALLY ONTO THE PANEL. 


Q 




PARABOLIC TYPE SOLUTIONS 

UNPROTECTED AREA 

ROD LENGTH IN METERS UNPROTECTED AREA IN % X ENHANCEMENT FACTOR 


.185 

1.1% 

2.9% 

.57 

.55% 

1.5% 

.74 

.28% 

.74% 


FIGURE 9 
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To illustrate the eone of protection concept for this configuration wo 
osg as an an example, a° 40°, 8 q 45°, and l o 7.3§m { *> 1/2 billboard 
width). The resulting Intersection Is a portion of an ellipse and Is shown on 
Figure 12. Even If those long (7.36m) lightning rods were placed every 
14.69m, a slgnlf leant fraction of the reetenna (6.7% or when weighted by 
enhancement factor 18%) Is unprotected (l.e. Is not Inside a cone of protee» 
tion). 

Furthermore, there are serious mechanical problems ass Jated with sup- 
porting those long (l.e., over 22 feet) lightning rods. We think these 
examples are sufficient to demonstrate that configurations employing fewer 
lightning reds at longer spacing decreases protection and creates structural 
problems that ultimately will Increase the total materials requirement. 

For example. If we were to Increase the length of the lightning rod In 
this configuration to the point that It could offer protection to the 
billboard In front of the one on which It Is mounted ( l.e. to the south), 
then with the appropriate phasing of rods between rows of billboards we could 
get total protection In the cone of protection context. The length of the 
rods would need to be 12m In order to provide this coverage. 

2.3 The Distributed Lightning Protection System 

The distributed lightning protection approach replaces the many lightning 
rods with a continuous horizontal conducting structure, as depicted In Figure 
13. The region of protection now becomes the volume beneath two planes whose 
Intersection Is the horizontal protecting structure. This protection tactic 
Is essentially the one employed by the power transmission companies. The 
angle between the protecting planes and vertical Is variable; 45° Is thought 
to be adequate but some designs use 30° for extra protection. This line is 
called the "static" by the power companies and this term Is used here for 
convenience. 

Figures 7 and 8 provide the correct geometric considerations for the 
distributed lightning protection If we Interpret the end point of the 
lightning rod to be the location of the static. We note that the figures 
apply anywhere along the reetenna, not just in the specific locations required 
by the lightning rod analysis. 

For consistent comparisons with the other lightning rod systems we will 
use a a 45 °. Since a < 45 0 for rectennas below 40° latitude, the top edge of 
the reetenna Is protected by the static for any value of L, the displacement 
distance. If we try to use a smaller, more conservative value for B, we will 
run into problems in protecting the top edge of the reetenna with any system 
tht does not cast a radio shadow on an active reetenna surface. The design 
constraint that we will use to specify L will be that the southward plane of 
protection intersect the reetenna surface at the base. Therefore, 

L ° 12.2m tan (45 0 - a). 

For a In the range 45 e to 30®, L has the range of values Dm to 3.3m. This 
simple analysis ignores the protecting capability of the Immediate southward 
row of the reetenna on the base of the row being considered. When these 
additional protective effects are considered we find that: 
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I o 6.1m (l - tan a) 

For a In the range 4S 0 to 30°, L now has the range 0m to 2.6m. 

Figure 13 gives the configuration of the distributed lightning protection 
system for a° 30®, whleh represents the most difficult situation to 
protect. In this situation the static Is displaced by 2.6 meters from the top 
edge of the reetenna; note that the 4S° planes of proteetlon provide total 
coverage of the reetenna. 

Wo wish to emphasize that the set of horizontal statics not only provide 
total protection In the sense that lightning flashes are expected to hit the 
statics Instead of the active reetenna surfaces but that this system also 
reduces the Induced voltages and currents In the reetenna. We estimate that 
Induced charges, currents, and potentials are reduced by 1/2 by the static 
protection system. 
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FIGURE 11 


PANEL SCALE PROTECTION COMPARED TO BILLBOARD 
SCALE PROTECTION SHOWN AS IN FIGURE 8 EXCEPT 
HERE ON A BILLBOARD. 


1*1 






FRONT 


FIGURE 12 

DISTRIBUTED LIGHTNING PROTECTION SYSTEM 


SIDE 
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FIGURE 13 


DISTRIBUTED LIGHTNING PROTECTION SYSTEM 
ILLUSTRATING FORWARD AND BACKWARD PRO- 
TECTION FOR SMALL INCLINATION ANGLES 
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II. SIMULATIONS OF LIGHTNING STRIKES TO THE SPS RECTENNA WITH AND WITHOUT 
PROTECT ION _ ' 


A series of experiments were performed In our electrostatic test chamber 
with a scale model of the SPS rectenna. The experiments consisted of exposing 
the model rectenna to a series of high voltage discharges produced with a 
Tesla coll . 

The strikes to the rectenna were photographed using time exposures In a 
darkened room. A wire from the upper plate conducted the discharge to the 
vicinity of the model rectenna and provided us with a limited control over the 
area of the strike. This allowed us to keep the strikes near the volume in 
focus by the camera. 

Different areas of the model rectenna were protected by different 
systems, and one area was unprotected. The following paracraphs describe 
samples of these experiments: 

1. The Unprotected Rectenna 

Most of the strikes were to the upper edge of the billboard because of 
the larger enhancement factor at that point. Several strikes to the billboard 
face occurred. 

In Figure 14, we see two strikes to the unprotected billboard section, 
one of which is to the billboard face. Notice that these strikes are perpend- 
icular to the face when near the face; we would anticipate this because the 
equipctential lines are nearly parallel to the face here. 

In Figure 14, we also see for comparison the' three lightning protection 
systems modeled. To the left is the billboard scale system; to the right is 
the panel scale system; and behind the flashes is the distributed lightning 
protection system. 

2. The Panel -Scale Protection System 

The next three figures are examples of strikes photographed on the 
section of the model rectenna that was protected by the panel -scale lightning 
protection system. 

In Figure 15, we see two strikes on the same billboard, both of which 
terminate on the panel-scale lightning rods. 

Figure 16 shows two strikes from a different view going to two different 
billboards. The panel-scale protection system here is seen to protect only 
the front billboard. Protection is probably greater for real lightning 
because in our experiments we artificially bring the "leader tip" very close 
to the billboard with the wire. 

Multiple strikes to the panel -scale protection system are seen in Figure 
18. One of the strikes goes directly to the billboard face, this type of 
failure will occur in nature, but with lower probability than illustrated 
here. 
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3. The Billboard-Scale Lightning Protection System. 

Two sets of experiments were made with the billboard-scale lightning 
protection system. The one Illustrated In Figure 19 corresponds to rods of 
length 7.35m. (A second series of strikes were made with rods cut to one-half 
of this length, but these were photographed In color and are not suitable for 
this report.) Figure 19 Illustrates the capability of these long rods to 
direct lightning to the desired point. 

In Figure 20, we have a side view of a billboard-scale protector taking a 
strike and protecting the billboard-face. Figure 21 Illustrates the "hole In 
the armor” of the bl 1 1 board- seal e lightning protection system. Two flashes 
strike the protection system, but a third strikes the billboards between two 
protectors, as predicted In Figure 12. With real lightning this Is less 
likely to happen, but It can and will occur. 

4. The Distributed Lightning Protection System. 

The displacement distance of the static from the billboard was scaled 
from 0.74m to make It correspond to the height of the panel-scale protection 
system. Fewer fallures-to-protect were observed with this system but they did 
occur. With real lightning, they would be even less likely to occur. 

In Figure 22, we see two strikes to two different billboards from the 
side view. Figure 23 shows two strikes to the same billboards, which were 
rovided with a distributed lightning protection system. One strike is to the 
terminal support rod at the billboard edge, which Is the preferred point of 
strike. The other strike goes to the horizontal static line between the 
terminal support rods. 
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III. GROUNDING CONSIDERATIONS FOR THE PRCPOSED LIGHTNING PROTECTION SYSTEM 


T he? thundercloud charges induce a large surface charge on the rectenna 
helow the cloud; during the stepped leader period even larger surface charges 
are induced on the region below the leader tip. Most of the current flowing 
auring the return strokes of the lightning flash must be distributed by the 
grounding system to connect with the induced surface charges. If adequate 
paths 'or these currents are not planned and provided, the lightning will make 
its own paths. Most of the induced surface charge will reside on the horizon- 
tal statics of the recommended distributed lightning protection system. The 
primary grounding system described here is to provide low impedance paths for 
the redistribution of the induced surface charges and the part of the light- 
ning charge that resides on the rectenna surface. 

1. Primary East-West Grounding 

It is absolutely necessary that the horizontal statics have a good low 
impedance connection at billboard edges. The static should appear to be a 
continuous very low impedance conductor in the east-west direction, as illus- 
trated in Figure 24. 

2. Primary North-South Grounding 

It is also necessary that the statics are mutually grounded in the north- 
south directions; there are two methods of achieving this: 

2.1 Periodic connections north-south at the level of the statics. If 
these north-south statics are aligned along the billboard edges, then 
there will be little power loss due to microwave shadows (See Figure 24.) 

2.2 Interconnect grounding in the north-south direction at the surface 
or sub-surface level (see figure 25) can also be used, but this approach 
creates a higher impedance to north-south currents on the static system. 

2.3 A surface level grounding network is required in addition to the 
primary static grounding network. The surface network must handle the 
redistribution of induced charges on the rectenna surfaces and power 
distribution systems and it provides a safe working environment at the 
surface level. Fast-west continuity with low impedance connections must 
be provided at the base of the rectenna support structures, and north- 
south continuity with low impedance connections as discussed in 2.2 and 
illustrated in Figure 25 must be provided. Figure 26 highlights the 
surface level grounding network. 

2.4. Interconnections between the primary and surface grounding networks 
should be provided by the vertical conductors located at every billboard 
upper corner; these are t he same structures on which are mounted the 
terminals and supports for the statics. The vertical interconnections 
are highlighted in Figure 27. 

2.5 The ultimate or final component of the grounding system is the tie- 
in to Earth ground. At regular intervals in the rectenna a deep earth 
qroundinq rod .must He cirjven intothe soil to make good contact with a 
conducting soil Tor earuh ground, ST 



f 


The organization of the earth grounding system should be along diagonals, 
as illustrated in Figure 28. Here we see that the placement of earth ground 
at every fourth billboard but on a diagonal produces a grid such that 
lightning striking the primary grounding network will never have to travel 
more than 30 meters along the east-west conductors before finding a ground, or 
32 meters along the north-south conductors (for a rectenna with a 40 5 inclina- 
tion angle). 
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GROUNDING 

RECTENNA LIGHTNING ROD SYSTEM 



FIGURE 25 
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Vertical interconnections between primary 

AND SURFACE NETWORKS AND ON EARTH GROUND 



FiniiRF ; >7 


33 



Placement of Eapth Grounds 



FIGURE 28 
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IV. MATERIALS AND SPECIFICATIONS FOR LIGHTNING PROTECTION 


It is premature to specify the final form for the materials for the 
lightning protection system. We think that the system should be integrated 
into the structural design of the rectenna itself; in this case many other 
considerations are necessary in addition to the capability to conduct 
lightning currents. The data displayed in Figure 2.9 (H. Baatz, Protection of 
Structures, in Lightning Vol . 2, ed. by R.H. Golde) i? useful for order-of- 
magnitude estimates of the lightning current requirements. 

Example: If the design permits a 100° C temperature rise in an ? aluminum 

member carrying 10° Amps for 10° seconds, we need approximately 3 mirr crossec- 
tional area of aluminum material in the conductor. pNote that the recommended 
crossections for building codes are larger 80 mm 4 ) indicating designs for 
lower temperature operation plus safety margins. 

The lightning conductor need not be solid. From a structural point of 
view a tubular or other extruded shape would be preferable. Such configura- 
tions are compatible also with the lightning protection recommendations. 
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Specific values of materials for wire 


Material 

Steel 

Copper 

Aluminium 

Density (g/cm -3 ) 

7-7 

8-92 

2-7 

Electrical resistance (12 mr ! m~*) 

0-17 

00178 

0029 

Heat (cal °C- 1 g- 1 ) 

0-115 

0093 

0023 

Melting point (°Q 

1,350 

1,083 

658 



. T , CmpCraturc ris = of conductors as fcnctoio of current square impulse per 
-oection square: Cu = copper. A! = aluminium, ST = steel, CRS = P 


cross-section square; 
rosion-resistant steel 


cor- 


Cross-scction for lightning conductors 


Installation 

components 

Material 

Cross-section 

(mm 3 ) 

Dimension 

Rod Strip 

(mm, radius) (mm x mm) 

Air termination 

Steel, galvanized 

50 (25)“ 

8 

20x2-5 

Rods up to 0-5 m 

Steel, stainless 

110 

12 

30x3-5 

long 





Down conductors 

Copper 

50(16)“ 

8 

20x2-5 

Conductors in 

Aluminium* 

80 (25)” 

10 

20x4 

ground 





1 

'Steel, galvanized 



0-5 mm 

Sheet metal J 

Copper 



0-3 mm 

1 

Aluminium, Zinc 



0-7 mm 

1 

Lead 



2.0 mm 


° Lowest cross-sections used in some countries. 
* Not for use below ground. 


FIGURE 29 
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v. ESTIMATE OF POWER LOSS FROM THE BEAM 


A rough maximum estimate of the power loss from the microwave beam due to 
the lightning protection devices can be obtained by assuming that the micro- 
wave shadow cast by the static lightning protection system is twice the 
crossectional area of the devices. We assume that the conductors are 2 cm 
wide of 1 mm thickness tubular material, providing 60 mrrr of crossectional 
area for conducting. The assumed shadow of these structures is approximately 
0.6% of the rectenna area (see Figure 30.), This is a maximum estimate of the 
loss. 
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VI. MICROWAVE DIODE FAILURES DUE TO INDUCED CURRENT TRANSIENTS 


The 25 W S GaAs diodes used in the design of the SPS reetenna have not 
been produced and no failure data is available for these devices. In order to 
obtain estimates of failure power of the diodes in the design, we used the 
specification data for the HP5082-2824 microwave diode and scaled the charac- 
teristics to 25 W using the "Wunsch relationship" described in the references 
below. We also obtained advice directly from Dr. D.C. Wunsch regarding the 
extrapolated power failure current. 

1. Defense Department Report D224-13042-1 EMP, Susceptibility of 

Semiconductor Components, dated September, 1974. 

2. Defense Department Report D224-10022-1 EMP, Electronic Analysis 

Handbook, dated May, 1973. 

3. Defense Department Report D224-10019-1 EMP, Electronic Design 

Handbook, dated April, 1973. 

Figure 31 shows the predicted failure power for 25 watt diodes, as a function 
of pulse width. 


39 



PULSE WIDTH - MICROSECONDS 


FIGURF 31 
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VII. COMPUTER SIMULATION OP ELECTROSTATIC FIELD AROUND AN SPS RECTENNA 


The electrostatic fields produced by the charges on the lightning channel 
Induce charges on the rectenna and on the lightning protection conductors. 
Changes in this electrostatic field require a redistribution of charge on the 
rectenna system; the resulting currents can cause diode failure even with a 
lightning grounding system in place. One output of the computer simulation of 
the electrostatic field around the SPS rectenna is an evaluation of the 
induced current on the rectenna with and without the recommended lightning 
protection equipment. 

An additional output from the computer simulation is the potential around 
the rectenna billboard enabling us to estimate the enhancement factors of the 
electric field due to the billboard shape. 


The algorithm used in the simulation computes an array of values for the 
potential around the middle of five infinitely long billboards. We assume 
here that the contribution to the local potential from billboards further away 
is ignorably small. The surface charge distribution on the billboards is 
simulated with ten infinitely long line charges evenly spaced along the bill- 
board. The value for the line charges is determined interactively with the 
computer to produce a zero potential contour that has the same shape as the 
billboard. Figure 32 illustrates this simulation. 


In order to compute the potential, we will need U(x,y), the electrostatic 
potential at a point (x,y) in free space, where the coordinate system is such 
that the line of electrical charges giving rise to the potential is located at 
the origin. If we call the y-coordinate the height h, then U(x,H) is the 
electrostatic potential at x and h of a line charge X (coulomb/meter) at a 
height d directly above the point x «* 0. There is also a contribution to U 
from the image charge. Thus, 

U(x,h) ■ in d ' | - g ] • 

d o x^ + (H + dr 

From this, the potential distribution around the rectenna may be calculated. 
Let U(l,h) be the potential at x * 1 and y = h due to a periodic system of 
line charges simulating the rectenna (see Figure 31.) We then have that 


uo.w At 


a_ „ n ,0 - ^ + fr- s h? m 

(1 - L[1 - 1] - X, f + (h + SX, f J ! 


2 ire 


J ' ' ’ ”‘j 

where the free-space value for the dielectric constant is assumed and where 


i = Billboard number, 
j = Line charge number on billboard i, 
s = Slope of billboard (= tan a), 

M = Number of line charges (= 10), 

N = Number of billboards (= 5). 
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SIMULATION OF SPS RECTENNA WITH LINE CHAR6ES 



FIGURE 32 



In the presence of a uniform electric field of 100,000 volts/meter 
(directeo upward), ten line charges have been selected to produce the array of 
values shown in Figure 33. Three potential contours have been sketched (zero, 
10,000 V, and 100,000 V) around the ten line charges on the billboard. The 
zero contour follows closely the position of the billboard surface, as 
required by the simulation algorithm. Note how closely spaced the contours 
are at the top edge of the billboard. Electric field enhancement factors of 
at least 6.5 exist in this region based upon our simulations. Higher 
resolution simulations would be required to refine the enhancement factor 
estimates. 

The values obtained for the 10 individual line charges found foi the 
solution shown in Figure 33 are (in p Coul./m): 

0.36, 0.465, 0.572, 0.679, 0.924, 1.02, 1.14, 1.78, 2.91, 4.14. 

We can convert these to a surface charge density by dividing each value 
by the billboard distance represented by the line charge. The first line 

charge serves approximately 3/2 n -)> the last line charge 

serves 1/2 ); and all others are associated with a length ). 

n 

Figure 34 is a plot of charge/unit area (p Coul./m ) on the billboard as a 
function of length (northward) along the billboard surface. 

When an additional line charge in placed at the position of the lightning 
static, and all of line charge values are adjusted to the new configuration, 
we find the simulated potential function around a protected billboard - Figure 
35. The placement of the static in this example is based upon the discussion 
in Section I. 2.3., with L =0.98m, corresponding to a = 40°. The charge/unit 
length for the static is 4.6 p Coul./m. The charge/unit lengths for the ten 
billboard line charges in (p Coul./m) are: 

0.315, 0.47, 0.51, 0.57, 0.87, 0.89, 0.90, 1.35, 1.78, 2.1. 

These line charges may be compared with the unprotected billboard charges 
corresponding to the solutions of Figure 35. The protected billboard charges 
approach approximately one-half of the corresponding unprotected charges. 

The line charges used to simulate the rectenna are normalized to a 
charge/unit area through division by the associated lengths, as previously 
described, to obtain the induced charge distribution on the protected rectenna 
bill board. 

2 

Figure 36 is a plot of charge/unit area in pCoul./m as a function of the 
distance (northward) along the billboard face. 
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ZERO POTENTIAL 


100,000 VOLT 

POTENTIAL 


10,000 VOLT 

POTENTIAL 


0 LOCATION OF LINE CHARGES SIMULATING BILLBOARD 
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FIGURE 33 
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VIII. COMPUTATION OF LIGHTNING ELECTRIC FIELDS 


In section VII, a rectenna was simulated in the presence of a uniform 
electric field of 100,000 Volts. The induced surface charges derived from the 
simulation are directly proportioned to the imposed electric field strength. 

In this section we describe a computer program that was written to derive 
values for the lightning-produced electric fields as a function of time and of 
distance from "ground zero" - the point of strike. We have run the program 
for a range of lightning parameters obtained from actual measurements reported 
in the literature. 

The program computes the contribution to the electric field from the 
thundercloud charge center participating in the cloud-to-ground flash, the 
charge on the lightning channel, and the images of these charges. All charges 
are allowed to vary with time in a manner consistent with observations [Ter- 
restial Environment (Cl imatic) Criteria Guide! ines for Use in Aerospace 
Vehicle Development, 1977 Revision ; Edited by John W. Kaufman, NASA Technical 
Memorandum ?8'ii8J. 

Figure 37 displays the relevant equations and configurations covering the 
leader phases of the computation. 

In Figure 38_the equations and conditions during the return stroke por- 
tion are shown. The program used in computing the fields is provided in the 
appendix. 

The material following Figure 38 provides the tabular and graphic data 
used in these computations for the return stroke phase. These data are 
contained in Figures (39-44) inclusive. 

The output of the computer program is a "blow-by-blow" history of the 
electrical field at a specified distance from ground zero as a function of 
tine. Figure 45 displays one section of the output from one of the computer 
runs. This corresponds to a worst-case situation, 10 meters away from the 
very-se\ere-model . The units of time are seconds( along the abscissa), and the 
units of the ordinate are kilovolts per meter. 

Table 8.4 in figure 46 provides a summary of the output for the various 
computer runs. Listed are the peak negative fields, the peak positive fields 
(when positive fields occur), and the £ and AT for the portion of the flash 
with the peak rate of change of electric field. 

These values are our input data to the computation of diode failure when 
used in conjunction with the induced surface charge results of the rectenna 
electrostatic simulations. 
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INITIAL SPECIFICATIONS 
TEMPORAL FUNCTIONS: 

X = Y 0 - Vj 
Q = Q 0 - P L (Y - X) 


(~5 km) 

Q 0 ( — 10 Coui) 

V L (— 10 5 m/s) 

t = o } Y = Y 0 
Q (~-5 Coiil) 

.P “ P L = Ql/Y-o 


SOLVE FOR E l (TjD) FOR t l WHERE 

t l = -V 7 \ 

X L (~50 METERS) 

FOR t>t l j E L (t,d) = E u (t l ,d) 



FIGURE 37 




P \ FOR Y >Y 0 MID P>-Pj_ 

Q = 5'IdT+ Q l J and g < — CQq - Q l ) 

SOLVE FOR E r (t',d) FOR t > t l OR t' ^ o 

TOTAL FIELD E t (t,d) = E L Ct l ,dO + E R (t,d) 

TERMINATE COMPUTATION WHEN Q - -CQ 0 - £ L ) 



v 2' QY' 

? 9 3/2 

Air £ q (D^ + Yq ) 


FIGURE 38 
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DETAILS OF A VERY SEVERE LIGHTNING MODEL (MODEL l) 


Stage 

Key Points 

Rate of Current Change 

Charge Passing j 

1. First Return 
Stroke Surge 

t= 0 i = 0 

t = 2 ps i = 200 IcA 

i 

t = 100 ps i = 7 IcA 

| Linear Rise - 100 IcA/ps 
| Linear Fall - 193 kA in 98 ps 

0. 2 C* 
~ 10.2 C 

2. First Stroke 
Intermediate 
Current 

t = 100 ps i = 7 kA 
t = 5 ms i = 1 kA 

| Linear Fall - 6 kA in 4. 9 ms 

19.6 C 

3. Continuing 
Current — 
First Phase 

t = 5 ms i = 1 kA 

t = 55 ms I = 400 A 

| Linear Fall - 600 A in 50 ms 

35.0 C 

: 

4. Continuing 
Current — 

Second Phase 

t = 55 ms i = 400 A 

t = 355 ms i = 400 A 

Steady Current 

120.0 C 

5. Second Return 
Stroke Surge 

t = 355 ms i = 400 A 
t = 355. 002 ms i = 100 IcA 
t = 355. 1 ms i = 3. 5 kA 

1 

]■ Linear Rise ~ 50 kA/ps 
j Linear Fall - 96. 5 kA in 98 ps 

~ 0.1 C 

~ 5.1 C 

6. Second Strcke 
Intermediate 
Current 

t= 355.1 ms i = 3. 5 IcA 

t = 360 ms i = 500 A 

| Linear Fall - 3 kA in 4. 9 ms 

9.8C 


* Coulomb (C) is the quantity of electricity transported in one second by a current of one ampere. 


FIGURE 39 


r--r 




























CURRENT INCREASING 


200 kA 


100 kA 



100 /is 


55 ms 


TIME INCREASING — 


355 ms J 355.1 ms 

355.002 ms 


DIAGRAMMATIC REPRESENTATION OF A VERY SEVERE LIGHTNING MODEL 
(MODEL l) (Note that the diagram is not to scale) 


FIGURE 40 
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DETAILS OF A 98 PERCENTILE PEAK CURRENT LIGHTNING MODEL (MODEL 2} 


Stage 

Key Points 

Rate of Current Change 

Charge Passing 

1. First Return 
Stroke Surge 

t= 0 1=0 

t = 5 jus i 100 kA 

t = 100 {is i 3. 5 kA 


Linear Rise - 20 kA/fis 
Linear Fall - 96. 5 kA in 95 ms 

0.3 C 
~ 4.9 C 

2. First Stroke 
Intermediate 
Current 

t = 100 ms i 3. 5 kA 

t = 5 ms i 500 A 


f Linear Fall - 3 kA In 4. 9 ms 

9.8C 

3. Continuing 
Current — 
First Phase 

t = 5 ms i 500 A 

t = 55 ms i 200 A 

1 

1 

i Linear Fall - 300 A in 50 ms 

17.5 C 

4. Continuing 
Current — 
Second Phase 

t= 55 ms i 200 A 

t= 355 ms i 200 A 

Steady Current 

60 C 

J 


FIGURE m 
























CURRENT INCREASING 



DIAGRAMMATIC REPRESENTATION OF A 98 PERCENTILE PEAK CURRENT 
LIGHTNING MODEL (MODEL 2) (Note that the diagram is not to scale.) 


FIGURE 42 
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DETAILS OP AN AVERAGE LIGHTNING MODEL (MODEL 3) 


Stage 

Key Points 

Rate of Current Change 

Charge Passing 

1. First Return 
Stroke Surge 

t= 0 i = 0 

t = 5 {is i = 20 kA 

t = 100 p.s i = 2 kA 


Linear Rise - 4 kA/fis 
Linear Fall - 18 kA in 95 pts 

0.1 C 
- 1.0 C 

2. First Stroke 
Intermediate 
Current 

t = 100 {is 1=2 lcA 

t = 5 ms i = 300 A 


i Linear Fall - 1. 7 kA in 4. 9 ms 

5.6 C 

3. Continuing 
Current — 
First Phase 

t = 5 ms i = 300 A 

t = 55 ms i = 100 A 


^ Linear Fall - 200 A in 50 ms 

10. G C 

4. Continuing 
Current -- 
Second Phase 

t = 55 ms i = 100 A 

t = 355 ms i = 100 A 

Steady Current 

30.0 C 


FIGURE 43 
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CURRENT INCREASING 



DIAGRAMMATIC REPRESENTATION OE AN AVERAGE LIGHTNING MODEL 
(MODEL 3) (Note that the diagram is not to scale.) 


FIGURE 44 
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FIGURE 45 
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TABLE 8, A 


VERY SEVERE MODEL 


98 PERCENTILE .MODEL 


AVERAGE MODEL 


Distance 

Peak 

Negative 

Peak 

Positive 

LY./L1 

Peak 

Peak 

Negative 

Peak 

Positive 

LE/&T 

Peak 

Peak 

Negative 

Peak 

Positive 

m/LT 

Peak 

10 m 

-8. 5X10 5 

2.8X10 5 

2.2X10 6 

-5. 95X10 5 

1.81X10 6 

6.46X1C* 

-5. 09X10 5 

1.39X10 6 

S. 68X10 5 



1.2X10" 5 



3.00X10 -6 



2. 59X10" 5 

50 m 

-5.7X10 5 

1. 7X10 5 

4.37X10 5 

-3.88X10 5 

1.04X10 5 

3. 59X10 5 

-3.1QX10 5 

6.1X10 4 

1.34X10 5 



2. 2X10" 5 



2. 5X10" 5 



2. 50X10" 5 

100 m 

-3.49X10 5 

2. 49X10 4 

2.15X10 5 

-2. 36X10 5 

N/A 

1.75X10 5 

-1.85X10 5 

N/A 

5.47X10 4 



2. 2X10" 5 


2.5 X10" 5 


3. 5X10" 5 

500 m 

-8.94X10“ 

N/A 

3. 79X10 4 

-6.15X10 4 

N/A 

2.96X10 4 

-5.12X10 4 

N/A 

N/A 


3.2 X10" 5 


4.5X10" 5 


1000 M 

-5.35X10 4 

N/A 

1.69X10 4 
4. 2X10" 5 

-2.61X10 4 

N/A 

N/A 

-3.29X10 4 

N/A 

N/A 


FIGURE 46 
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IX. COMPUTATIONS OF DIODE FAILURE 


We are now to the point of having generated all of the data that are 
required to evaluate the conditions under which the microwave rectifier diodes 
will fail due to induced currents from nearby lightning flashes. For a 
given /£ and AT (from Table 8.4) we obtain from Figure 31 the power required 
for diode failure and from Figure 32 the induced charge/unit area on the 
rectenna surface. We assume that a diode designed to operate at 67 V will 
have a breakdown voltage of about 100 Volts. 

The surface area of the rectenna that has an induced surface charge of 
the size sufficient to cause diode failure is then computed from comparison 
with areas of the rectenna served by individual diodes and by series strings 
of diodes. Sample computations follow. 

SAMPLE COMPUTATION OF DIODE FAILURE 
(98TH PERCENTILE - 10 METER - NO PROTECTION) 

1. 98 percentile model - 10 meters: AT = 3 x 10" ®and t£ = 6.46 x 10®. 

2. Expected diode failure power from Figure 30: 250 Watts. . 

3. Energy dissipated in the diode: 250 Watts x 3 x l" b s = 7.5 x 10"^ 

Joules. r 

4. Charge transferred across 100 Volts diode breakdown voltage « 7.5 x 10 _t> 

Coulombs. 

5. From t£ ip step 1 and figure 37, the induced charge/unit area 

= 3 xlO ' 6 c/ iii x 6.46 = 19.38 x 10 “ 6 c/nf. 


6 . From steps 4 and 5, the rectenna area with surface charae equivalent to 

the charge required to cause diode failure is: 0.39 nr. 

7. Area served by diodes: rectenna center, 


25 watts „ ,, 2 , 25 watts „ r 2 

—o- s 0.11 m ; rectenna edge, 5 — = 2.5 m . 

230 w/ nr 10 w/nf 

8 . Compare 6 with 7: single diode configuration near rectenna center is 

safe. Single diode configuration near rectenna edge is vulnerable. 

9. However , the diodes are to be put in series (597 to a string) hence the 

diodes near the bottom must carry all of the induced current to the 
entire string. For these bottom- string diodes the area served with 
respect to the induced charge is: rectenna center, 60 nr; rectenna edge, 

1400 m . 


10. To protect against the 98 percentile flash within 10 meters of ground 

zero would require fast surge protection diodes (back to back zeners) on 
all diodes in the rectenna. This extent of protection may not be cost 
effective; however the considerations in Section X indicate that simpler 
protection arrangements will probably be effective near the rectenna 
center. 


FAILURES PRODUCED BY THE AVERAGE LIGHTNING FLASH 


The situation considered here is the extent of the protection required 
for an "average" lightning flash if we are willing to accept losses from the 
extreme cases. 
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The computation sequence follows the same procedure described immediately 
above. Here we use data for the average flash from Table 8.4 at a 10 m 
distance from ground zero. 


SAMPLE COMPUTATION OF DIODE FAILURE 
(AVERAGE FLASH, 10 M, WITH "STATIC" PROTECTION) 


1. From Table 8.4: £ = 5.68 x 10 5 v/m; AT = 2.59 x 10" 5 s. 

2. From Figure 6.1: 80 watts. o 

3. 80 w x 2.59 x 10 s « 2 x 10" ^ Joules. 

4. 2 x 10" b coulombs. r , 9 

5. From 1 and Figure 38: 1.5 x ? 10~ 0 x 5.68 *» 8.52 x 10" coul/m. 

6. From 4 and 5: Area = 2.35 m. 

7. Since the rec|enna area served by individual diodes even on the 

edge < 2.5 m , the individual diodes are self- protecting and able to take 
an "average" lightning flash. 

8. However , when arranged in a series stack of 597, the diodes at the 
bottom of the stack must conduct the induced currents for the whole 

stack. The diodes cannot safely carry these currents. 
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X. LIGHTNING PROTECTION FOR SERIES DIODE STRINGS 


As demonstrated In Section IX, the connection of microwave rectifier 
diodes in series requires special lightning protection considerations. We 
cannot make specific recommendations for these protection devices at this time 
because the rectenna current design is not advanced to the point that allows 
such detailed analysis. Rockwell International has provided us with an equiv- 
alent circuit for the rectenna; a slightly modified form of that circuit is 
shown in Figure 46. We have assumed that the series connections are to be 
made at the points indicated by the large spots and that the output filler 
operates around 30 Hz. A series string of rectenna elements of this design 
can be protected with a variety of methods. One cost-effective means is a 
spark gap arrangement incorporated in the diode feedthroughs, or the output 
filter inductors, or on the billboard configuration itself. 
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RECTENHA EQUIVALENT CIRCUIT AT 2,45 GHz 



FIGURE 46 
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XI. CLOUD -TO-GROUND LIGHTNING DISTRIBUTION IN THE UNITED STATES 


In order to have a working estimate of the hazard presented by lightning 
to rectennas, we need to know the cloud-to-ground lightning flash density for 
various possible rectenna sites in the United States. The cloud-to-ground 
lightning flash density (in #/knr for example) is not a parameter that is 
measured as a climatological variable. We have found it necessary to use the 
number- of- thunderstorm days as a proxy variable because it is available as a 
climatological variable. Figure 47 gives contours of annual number-of- 
thunderstorm days. 

XI. 1. Pierce Conversion Formula 

Several attempts have been made to derive a conversion formula to convert 
thunderstorm days into the flash density by using lightning flash counters in 
research areas for correlation with the count of thunderstorm days. The best 
of the various conversion formulas is that due to E.T. Pierce ("A Relationship 
Between Thunderstorm Days and Lightning Flash Density," Trans. AGU, 49 , 686, 
1967.) The Pierce formula (as does most others) has a quadratic term, which 
reflects the relationship between frequencies of local storms and storm 
intensity. In addition, the formula utilizes the monthly thunderstorm days as 
opposed to the annual average in order to incorporate seasonal effects in the 
conversion formula. 


This formula is 


? ? 4 

= aT M + a*V 


'M 


■M 


where: T M = monthly number of thunderstorm days and a, is the monthly ground 

flash density (#knr/Mt.) The parameter a is, p 

a = 3 x 10" 2 

If a is the annual ground flash density (// knfVyr.), then 

c = V* o . 

(M M 


XI. 2. Climatological Data -- Number of Thunderstorm Days 


The inputs needed to compute the U.S. Distribution of ground lightning 
flash density are: (1) The monthly number of thunderstorm days for all U.S. 
stations recording these observations, (2) the coordinates of the observing 
sites, and (3) the computer software to compute the density and display the 
results geographical ly. 

Items 1 and 2 were obtained from "Local Climatological Data - Annual 
Summaries for 1977" published by The National Oceanic and Atmospheric Adminis- 
tration on magnetic tape. The geographic plotting software of Item 3 was 
obtained from The National Technical Information Service, and the computer 
programming was done by J.L. Bohannon at Rice. 
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A detailed list of flash density for all of the stations used is provided 
in the Appendix. 

Note the hot spots on the contours in Figure 48 that result when stations 
are located near geographic features that promote local thunderstorms. There 
are probably other similar hot spots in the U.S. that do not show up on this 
display because of the absence of an observing station nearby. 
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FIGURE 43 
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APPENDICES 


Computer programs developed under this contract. 

All programs are in FORTRAN H, unless otherwise specified. All of the 
programs were run on an IBM 370/155 and/or an Itel AS/6 computer. 
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A 


Computer Program PANEL: 

A Computer Model of the SPS Plasma Interaction 


The following pages are the listing of the program "PANEL," written to 
model the interaction of a high voltage solar array with an ambient Maxwellian 
plasma. The program was originally written by Dr. Lee W. Parker and was 
modified for application to the SPS problem by David L. Cooke. 
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PRINTERS 

A 

END 

JOB 

2 AO A 

BIN 

7 

COOKE 

437075SCP5 

4:02:00 

PH 

27 

HAR 

79 

PANEL 


PRINTERS 

A 

ENO 

JOB 

2 A 0 A 

BIN 

7 

COOKE 

A37075SCP5 

A: 02:00 

PH 

27 

HAR 

79 

PANEL 

J 

1 PRINTERS 

A 

END 

JOB 

2A 0 A 

BIN 

7 

COOKE 

A37075SCP5 

A: 02: 00 

PH 

27 

HAR 

79 
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PRINT ER3 

A 

END 
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END 

JDB 

2 A OA 

BIN 

7 

COOKE 

A37075SCPS 

A: 02 : 00 

PH 

27 

HAR 

79 

PANEL 

V 

^ ft 

END 
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PH 

27 
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79 
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END 
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7 
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PH 

27 
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79 

PANEL 


NT Ek 


END 
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PH 

27 
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79 
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27 
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79 
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3r 7-rj 
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LEVEL 21.8 C JUN 74 ) 


QS/360 FORTRAN H 


DATE 


I 


COMPILER OPTIONS - NAHE= HAI N, OP T=0 2 , L IME£NT=60. S I ZE=Q OOOK, 

SDURC E, EBCDIC .NOLIST , NOD ECK, LOAD, KA P ,NOEOIT, HQID , NQXREF 
C 

C SOLAR PANEL PROBLEM 


ISN 

0002 


ISN 

0003 


ISN 

0004 


ISN 

0005 


ISN 

0096 


ISN 

0007 
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00 0 8 
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0009 


I 5N 
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QG 1 2 
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0013 
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I 5 N 

0015 

2 

ISN 

0016 
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C 

£ 
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rsw 
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C 
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0032 
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003 3 
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0034 
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0035 


ISN 

0036 


T 5 N 

003 7 


I SN 

00 3 H 


I SN 

00 39 
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C 

C 

rsN 

0040 


I S3. 1 

0942 


I SN 

004 3 


I SN 

004 4 


1 SN 

00 4 5 


ISN 

0 0 4 6 


1 SN 

0047 



C0M u QN/CP/NPRINT,NPTS,KA*KB,HE*KHAX,XPT*YPT»ALl»8El*EVfSHACH* 
1 T VOLT St CUR, XHE TER 

COMHON/BK/IIH, IIP, JJH, JJP, KK,NTOT,IV, JV, II, J J,K* N. VP C30) , 
1XYZC2G8G »3)»VVC3Q»20,1G) , XP C 30) , XKC 1 0} , YP C 203 , YKC1Q) , ZLC 1 0) » 
2XXC40) ,YYC30) , ILX , IUX * .<UK , HOC , Ha 0 * VRF, NFP S * SKPRFL, SKPLST 
COMM')N/FLD/XC2089,2),CO£FC208Q,7),INDXC298Q,6),SKPC0 
COMMON /CO/P VOLT 5 , XM ACH , OE NS T , NN , P AR7CLC2 ) , PARTI (23, PART2C2) 
COMMON/ IN TER / INT * 1 1 A * J JA , KKA , 1 GOUT, JGOEJT , KGOUT»XA*YA , ZA , 
1XI(30),YJC20),ZKC10) 

DIMENSION DATE C20> 

DIMENSION VrC C4), IFC 4 ), JFC4),KFC4 3 

INTEGER S<PRFL,SKPLST,S<PCO 

NFC IX » JX ,XX)=JX* II*CJX-l)*-II*JJ*OCX-i) 

L=5 
M= 6 

READCL,9999, END=99) DATE 
FORMA T C2GA4) 

WRITE CM, 9998) DATE 

FORMAT C42H1S0LAR PANEL ELECTRIC FIELD AND CURRENTS. ,20A4) 
READ GEOMETRIC PARAMETERS 
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READ PANEL POTENTIALS 
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N = NFCT ,J,1) 
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CONTINUE 


CONSTRUCT REFLECTORS 
IFCSKPRFL.EO.DGO TO 
DO 169 I = ILX, IUX 
DO 160 K = XLK.jCUK 
J A ^ HBC-X 
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C 

163 

118 


119 

117 

165 


170 

220 


C 

111 

113 


JW = K+HBD 
NW = NFC I , JW , K 3 
XCNW.23 = 1 
XCNW.13 = VRF 
CONTINUE 

WRITECM, 2313 VRF 
FORMAT C//1X , * REFLECTOR 


POTENTIAL = * *1PE15 .5) 


FIXED POTENTIALS 
TO 220 


FIXED POTENTIALS*/ 
i 3 X » ‘ K 1 33 


READ ADDITIONAL 
IF CNFPS.LE.33G0 
WRITECM, 1183 
FORMAT C//' ADDITIONAL 
1 4 C 6X , * POT * 

DO 170 NOC = 1,NFPS,4 
READ CL, 119 3CVFC(I3,IFCI3 , JFCI3 ,KFCI3, 1=1 , 43 
F0RMATC4CE8. 0,31433 

WRITECM, 1173CVFCCI3.IFCI3, JFCI3.KFCI3. 1=1,4 
FOR MAT C/4 C3X.1PE1 0.2 ,314 33 
00 170 1=1,4 

HN = NFCIFC 13 , JFCI3.KFCI33 

XCNN, I 3=VFCCI3 

XCNN, 2 3=1 

CONTINUE 

CONTINUE 

I VP=TV +1 

JVP=JV+1 

WRITECM, 113 3 IIP, IIM, JJP, JJH,KK, IV, JV 
WRITECM, 2233 Cl , XP Cl 3 , 1 = 1 , I V3 

CI. XPCI3 ,I=IVP,IIP3 
CI«XMCI3,I=1*IIM3 

CJ, YPCJ3, J=l, JV3 
C J, YPCJ3, J=.JVP, JJP3 

CJ, YMCJ3 , J=l, JJM3 

CK, ZZCK3 ,K=1,KK3 
CXPCI3 , 1=1 , IV3 


WRITECM, 2243 
WRITE CM, 225 3 
WRITECM, 2263 
WRITECM, 2273 
WRITECM, 2283 
WRITECM, 2293 
WRITECM, 2303 


O 7, 

|* 

^ 03 


WRITECM,2413CVPCI3,I=1,IV3 


1 

2 

3 

4 

5 

6 


FORMAT 

FORMAT 


ISN 

0084 

116 

FORMAT 

TSN 

0085 

222 

FORMAT 

ISN 

0086 

223 

FORMAT 

ISN 

0087 

224 

FORMAT 

ISN 

0088 

225 

FORMAT 

I SN 

0089 ■ 

226 

FORMAT 

ISN 

0090 

227 

FORMAT 

ISN 

009 1 

220 

FORMAT 

I SN 

0092 

229 

FORMAT 

ISN 

909 3 

230 

FORMAT 




1 15X.3 

ISN 

0094 

240 

FORMAT 

ISN 

0095 

241 

FORMAT 


Cl 61 53 

C//1X, I3.18H 
IX, 13, 18 H 
IX, 13,1811 
IX, 13, 1811 
IX, 13, 2511 
IX, 13,3311 
IX, 13 , 33H 
C8E10.03 
C16E5.03 

C//1X, 27HX-VALUES 
C//1X,35HX-VALUES 
C//1X, 18I1X-VALUES 
C//1X.27HY-VALUES 
C//1X,35HY-VALUES 
C//1X, 10HY-VALUES 
C//1X,37HZ-VALUES 
C/ ///IX, 25 M ARRAY 


POSITIVE 
NEGATIVE 
POSITIVE 
NEGATIVE 
Z -VALUES 
POSITIVE 
POSITIVE 


X-VALUES/ 

X-VALUES/ 

Y-VALUES/ 

Y-VALUES/ 

CPOSITIVE ' 0NLY3/ 

X-VALUES DEFINING PANEL/ 
Y-VALUES DEFINING PANEL3 


MX =,3X,C8CF8. 4, 4X3/20X33 
C/1X,2HYC, 12,211 3 = , F8. 4 „6X 


POSITIVE ON PANEL=/CI3,lPEl5-433 
POSITIVE OUTSIDE OF PANEL=/CI 3, 1PE15.433 
NEGATIVE=/CI3,1PE15.433 
POSITIVE ON PANEL=/CI3,1PE15.4)3 
POSITIVE OUTSIDE OF PANEL=/CI3, 1PE15.433 
NEGATIVE = /CI3, I PEI 5 -4 33 

CPOSITIVE ONLY 3 ABOVE P ANEL=/CI3, 1PE15.433 
OF PANEL POTENriALS// 


C8X, ' ALL Y, • ,5X, C8C1PE 12. 43/20X33 


C8CIPE12- 43/20X33 


PAGE 003 


I , 1 5 1 5 ) 

Y1X,3',HNPRINT,NPTS,MA,MB,ME,KMAX ,PR0BN0 =,616,110/ 

Nil MB £11 =, F9.l t 9X, 1311 TEMPERATURE =, F9.1, 6H VOLTS, 9X, 
F'J.l, 7U PER CC, 9X, 6MMASS =, F9. 0,1111 ELECTRONS / 
■ T I SCALE = « F9 . 1 , 30M METERS = X-DIMENSION OF PANEL! 

? JllSINGLE SPACE POTNT. X =,F10-5,5X, 3HY =,F10.5! 

3AM SINGLE ENERGY ( MONOENERGET IC 1 . E =,F10.5, 6H VOLTS! 
22TSTNGLE TRAJECTORY. X =,F10.5,5X, 3HY =,F10.5/ 

LI ALPHA = , F20 . 8 , 8H DEGREES/ 

■LL META =, F20. 8 , 8M DEGREES/ 

R GY =, F20.5 , 6M VOLTS! 

32MRANOOM THERMAL CURRENT DENSITY =,1PE13.'A, 

FI SOUARE METER, F0R.2A5! 

II .INTERFACE X - VALUES / C I 3 , 1 PE 15 . A !! 

1 ll MINT ERF ACE Y -V ALUES /C T 3 » 1 P El 5 . A 1! 

1 Ull INI ERF ACE L- V ALUES / Cl 3 , IP E 1 5 . A TJ 

f " 22M -- CURRENTS AND POWER!! 


r i 


\ ' 

r uO TO A 20 


*'< I 
>' ! 


l.i 


b 


0 TO A 00 


' . J 
) 


1 \ 0 )- Yfl(i) - -5 

/Mr 1 *)( 1 Q 

7^0 


1 

i 


x I-lMXXCIl! 


/ 

rj-ll«YYCJ!l 


1 


y 



ISN 

0141 


ZKC13=ZZC13 

ISN 

01 42 


ZXCKKA3=ZZ(XX3 

tsn 

0143 


DO 560 K=2,XK 
ZKCK)=.5*CZZCK-13+ZZCK33 

1 SN 

0144 

560 

1 SN 

0145 


WRITECH, 5611 C I , XI Cl 3 , 1=1 , IIA3 

I SN 

0146 


WRITECM.5623 C J, YJCJ3 , J=1 , JJA3 

ISN 

0147 

r 

WRITE CM* 5633 CK , ZX CX 3 , K=1 , KKA3 

TSN 

014Q 

u 

DO 600 N=l*NTOr 

tsn 

0149 


CALL FINDCIFIND, JFIND .KFIND3 

TSN 

0150 


XYZCN, 13=XXCIFIND3 

ISN 

0151 


XYZCN, 23=YYCJFIND3 

ISN 

0152 


XYZCN, 33=ZZCKFIND3 

ISN 

0153 

600 

r 

CONTINUE 

ISN 

0154 

U 

IFCSXPLST.EQ.13 GO TO 660 

ISN 

0156 


NFPP = C NTQT/3Q03+1 

ISN 

0157 


DO 650 IP=1 , NFPP 
WRI TECH, 9OO03 

ISN 

0158 


ISN 

0159 

9000 

F0RMATClHl/6X,lHN,3X,4HXCN3,2X,4HYCN3i2X * 4HZCN3//3 

ISN 

0160 


CALL L ISTC2, IP3 

ISN 

0161 

650 

CONTINUE 

ISN 

0162 

660 

r 

CONTINUE 

ISN 

0163 

L 

00 700 J=1,JJ 

ISN 

0164 


DO 700 1=1,11 

ISN 

016 5 


N = NFCI.Jil3 

ISN 

0166 


VVCT.J.13 = X CN* 13 

ISN 

01 67 

700 

r 

CONTINUE 

ISN 

0160 

L 

K=1 

ISN 

0169 


WRITECH, 8000 3 K, ZZCK 3 , CXX CI3 , 1=1 » 113 

ISN 

0170 


OU 750 J=1,JJ 

ISN 

0171 


WR I TECH, 24 03 J, YYCJ3 , C VVC I , J ,K3 , 1=1 , 113 

ISN 

0172 

750 

r 

CONTINUE 

ISN 

0173 

lr 

r 

CALL FIELD 

TSN 

0174 

b 

00 800 K=1 , XX 

ISN 

0175 


no ooo . 1 = 1 , jj 

ISN 

0116 


DO 800 1=1,11 

ISN 

0177 


N=MFCI, J,K3 

ISN 

0178 


VVCI,J,K.3 = X(N, 13 

ISN 

0179 

800 

r 

CONTINUE 

ISN 

0180 

L 

DO 900 K=L,KK 

ISN 

0181 


WRITECH, 00003 K, ZZCK3 , C XXCI 3 , 1 = 1 , 113 

£ SN 

0182 

8000 

FOR HAT C 26M1ARRAY OF POTENTIALS AT ZC, 12 , 2JO=, F8. 4// 

1 1 5 X , 3UX =,3X, C8CF0. 4, 4X3/20X33 

ISN 

018 3 


DO 850 J=1,JJ 

ISN 

0104 


WRITECH, 24 03 J , YY C J 3 , C V VC I , J, K3 , 1=1 , 1 13 

1 SN 

0185 

850 

CONTINUE 

ISN 

0186 

900 

r 

CONTINUE 

ISN 

0187 

U 

NP ROB= 0 

ISN 

0188 

1000 

REAOCL ,333,END=993 M PRINT, NPTS, HA ,HD , HE ,KHAX ,HORE 
REA DC L ,1163 SHACH.TVOLTS ,DENCC,XHASS .XHETER 

ISN 

0109 

1001 

ISN 

0190 


NPROB = NPROBH 


1 SN 0191 
I5N 0192 

1SH 0193 
TSM 019 5 
ISM 0197 
I S II 0199 
I5N 0201 
I 911 0203 
I S 11 0205 
ISM 0207 
TSM 0209 
ISM 0211 
TSM 0213 
ISM 0214 
TSM 0215 
TSM 0217 
ISM 02 1 0 
ISM 0219 


WR1TECH,9991 

WRITE CM, 44 4 9 NPRINT, HPTS * HA » KB » HE , XHAX >NPROB »SHACH »TV0LTS ,OEMCC. 
1 XHASS ,XMETFR 

IFCNPTS.EQ.O.Ok.HE.EQ.O.QR.hA.EQ-09 RE AD CL , 2229XPT , YPT , AL1 , OEl , EV 
IFCMP1 S.E0.09 URITECH, 445) XP7*YPT 
IFCKE.E0.Q9 WRITEC’1,4 4 69 EV 

IFCMA.EQ.Ol WRITECM.4479 XPT, YPT ,AL 1 , BE1 , EV 
IF CMA.GT.O. AND. XHASS .LE.0.3 STOP 

IFCHA. GT.09 CUR=2.63£-8*DEMCC*SQRTCABSCTV0LTS)/XHASS) 

IFCTVOLTS.GT. 0.) PAR TCL C 1 )=PART 1 C 19 

IFCTV0LTS.GT.0.9 P AR TCLC 2 9=PAR T1 C29 

IFCTVOLTS.LT. 0.3 P AR TCL C 1 9=P ART 2 Cl 3 

IFCTVOLTS.LT. 0.9 PARTCLC2 9=PART2C29 

WRI TECH, 4489 CUR,PARTCL 

CALL POWER 

IF CHORE. GT. 09 GO TO 1000 
GO TO 100 
99 STOP 
EMO 



i. * t. 1_ 


CxtU V JUH [H j 


Uo/ibU MJKiKAN H 


UAfE 


ISN 0002 
ISN 0003 


ISN 

ISN 


0004 

0005 


COMPILER OPTIONS - NAME = M AIN, 0PT=02 , LINECNT=6Q , SIZ E=0000K , 

SOURCE, EBCDIC, NOLIST, NODECK, LOAD, HAP , NOEOIT, NO ID, NOXREF 
SUBROUTINE ORBIT 

STEP ACROSS 3-D BOX ASSUMING CONSTANT POTENTIAL WITHIN BOX 

COHNON/BK/IIH, IIP , JJH, JJP ,KK,NT0T,IV, JV, 1 1 . J J, H, N , VP C30) , 
1XYZ(2O8O,3),VVC3O,2O,1O),XPC3O),XHC1O),YPC20).YMC1O),ZZC1O), 
2XXC40) » YYC30 3 , ILX , IUX , KUK , HBC , HB D , VRF, NFP S , SKPRFL , SKPLST 
COMMON/ORB/XDOT, YOOT,ZDOT, XI, X2, Y1,Y2, 21,22, X^Yil, PHI, NT IHE, SAVE 
DIMENSION TIMEC6) ,UC3) , UDOT C3 3 , B C2, 3) 


ISN 0006 
ISN 0007 
ISN 000 B 

ISN 0009 
ISN 0011 
ISN 0013 


ISN 0014 
ISN 0015 
ISN 0016 


TSN 

ISN 

ISN 

ISN 
ISN 
ISN 
ISN 
I SN 
ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 


ISN 
ISN 
I SN 
ISN 
TSN 
ISN 


0017 

001U 

0019 

0020 
0021 
0022 

0023 

0024 

0025 

0026 
0027 
0029- 
0030 
003 1 

0032 

0033 

0034 

0036 

0037 


0038 

0039 

0040 
0042 
0044 
0046 


C 

c 

c 


T00H=3.3333E+33- 
ROUND = l.E-11- 
T R-GUN U — =— irrfr-Yj 


^JT F./ cT S ) 


999 

C 

C 


IFCXOOT-EQ.O.-ANO. YDOT.EQ. O-.AND.ZDOT. EQ.G. 3 WRITECH.999) 
IFCXDOT.EQ.O.. AND. YDOT. EQ.O.. AND. ZDOr.EQ.O .) RETURN 
FORMAT (1X,3BHSPEED=0 - HENCE PARTICLE DOES NOT HOVE) 


UC1)=X 
U C 2 3 = Y 
UC 3)=Z 

UDOTC 1 )=XDQT 
UDOT C 2 ) = YD0I 
UDOTC 3)=ZDQT 

B C 1 , 1 3 —X 1 
BC2,1)=X2 
B C 1 , 2 )=Y 1 
BC2,2)=Y2 
B(1,3)=Z1 
BC2, 33=72 


j£jL. 


1 

l f 

5 

0 




Y-l 


J 


o c 
^ 3? 
tj S 

o 5 

*D f— 
O TJ 

MS 

I" 

"S 03 


01 


on 101 N 2=1 , 3 

IFCUDnTCN2).EO.O.) GO TO 101 
DO 100 HI =1,2 
N R = N 1 + _2*XN2-1) 

T I M E CN R 3 = r 0 On 

'TT=C (I C N1 ,N2) J- UCN233/UDOT-CN2) 

SS = UC.N2) n- JJD0TCN2)*TT 

IFC SS 7'GC'^B‘T 1,N2).AND.SS.LE.BC2,N23) 

CONTINUE 

CONTINUE 


100 
1 
C 

C FIND SHORTEST SIGNIFICANT TIME 
C 


TIHECNR)=TT 



TIM IN=TOOH 
DO 200 NR=1 , 
IF CTIMECNR) . 
IFCTIHECNR). 
IFCTIMECNR). 
CONTINUE 



CrX' /{ 

tyvvJJ '■u s)C\j 


S5- 




& 


k 


200 
C 

C ADVANCE TO APPROPRIATE END-POINT 
C 


EO.TOOH) GO TO 200 ‘ 

GT-R0UND.AND.TIMECNR3.LT.TIHIN) .NTIHEpNR 

GT. ROUND. AND. TIMECNRO.LT.TIMIN) C TIHIN=TIHE CNR) 


? 


1 


ISN 0047 
ISN 0048 
ISN 0049 

ISN 0050 
ISN 0051 
ISN 0052 

ISN 0053 
ISN 0055 
ISN 0057 
ISN 0059 
ISN 0061 
ISN 0063 

ISN 0065 
ISN 0066 
ISN 0067 

ISN 0068 
ISN 0070 
ISN 0072 
ISN 0074 
ISN 0076 
ISN 0078 

ISN 0080 
ISN 0001 


C 

C 

c 

c 

c 


X=X + X00T*TIMIN 
*Y = Y + YD0T*TIMIN 
2 = 7. f ZDOT*TINlN 

XS A V=X 
YS A V = Y 
ZSAV=Z 


IFCNTIHE.EQ.l3 
IFCNTIME-EQ. 23 
IFCNTIME.EQ.33 
IFCNTINE.EQ.43 
IFCNTIME-EQ. 53 
ie.CNTIHE.EQ. 63 


Y=Y 1 
Y=Y2 
Z=Z1 
Z=Z2 


DX=X-XSAV 

DY=Y~YSAV 

DZ=Z-ZSAV 


Z' 


IFCCNTIHE.EQ.1.0R.HTIHE.EQ.23 
lFCCNTIHE.EQ.3-0R-NrlHE.EQ.43 
IFCCNTIHE.EQ.5.0R.NTIHE.EQ.63 
IFCNTIHE.EQ.-13 SAVE=XSAV 
IFCNTIME.EQ.-23 SAVE=YSAV 
IFCNTIHE.E0.-33 SAVE=ZSAV 


RETURN 

END 


* 


i 


i: -fc (X/ 




\ 


AND ABS CDX3 . 6T.TRQUND3, NTIHE=-1 
AND.. ABSCDY3.GT .TRQUNQ3.' NTIHE=-2 
AND.. ABSCDZ3.GT .TR0UND3 NTlHE=-3 

V . 






LEVEL ^l.U C JUN 7 4 3 


U5/36D FUKIKAN f! 


Ufl j c 


TSN 0002 


COMPILER OPTIONS - NAME = MAIN, 0PT=02,LXNECNT=60, SI ZE =0000 K, 

SOUR CE, E BCDIC, NOLIST, NOD ECX, LOAD, HAP, NOEDIT, N0ID,N0XREF 
SUBROUTINE DEN 


ISN 0003 


ISN 0004 


ISN 

ISN 

ISN 

TSN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

TSN 

ISN 

ISN 

TSN 

ISN 


ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 


ISN 

ISN 

ISN 

ISN 

ISN 

ISN 


0005 

0006 

0007 

0008 

0009 

0010 
0011 
0012 
0013 

0015 

0016 
0018 
0020 
0021 
0022 

0023 

0024 


0025 

0027 

0028 

0029 

0030 

0031 

0032 
0034 


0035 

0037 

0030 

0039 

0040 
004 2 


ISN 0043 


ISN 

ISN 

TSN 

ISN 


0044 
004 5 
0046 
00 4 8 


C 

C 

c 


ROUTINE FOR EVALUATING CURRENT-DENSITY INTEGRALS OVER VELOCITY SPACE 


CQHMON/BK/IIM,IIP,JJH,JJP,KK*NTQT, IV, JV, II , J J,H,N,VPC303 , 
1XYZC2080,33, VV(30, 20 , 10) , XP (30) , XMC 1 03 , YP ( 20 3 , YHC10) . ZZC103 , 

2XX C4Q 3,YYC3u),ILX,IUX,KUK,MBC»HOO,VRF, NM* *» . SKl'RFL, SKPLST 
CDHHON/CP/NPRINT,NPIS, MAjMB, HE , X.MAX , XPT , YP 1 ,AL1 , BE1 ,EV ,SHACH, 

1 TVOLTSjCURjXHETER 

COMMON /CO/P VOLTS ,X M ACH , 0 ENS T , NN f P ARTCL C23 , PART1C23 * PART2C 23 
COMMON/ORB/XDOT , YOOT , ZOOT , XI, X2, Y 1, Y2, Z1 , 22 , X, Y , Z, PHI, NTIME , S AVE 
COMMON/ INTER/ IN T ,IIA, JJA, KKA, IGOiUT, JGOUT ,KGOUT,XA, YA, ZA, » 

1XIC30) ,YJ(203,ZKC103 
DIMENSION AC23,END1C23 , END2C23 , FATEC 23 
OATA ENDl/41IABS0,4HRBED/,EN02/4KESCA,4HPES / 

XS A VE=XPT 
YSA VE= YPT 
TE MP= ABSCTVQLTS3 
IFCIEHP- LE . 0 - 3 WRITECH,9993 TEMP 
999 FORMAT C////1X, 3BHTR0UBLE - NEGATIVE OR ZERO TEMPERATURE 3 
IFCrEHP.LE.0.3 RETURN 
IFCMA-EQ.O.OR-HE.EO. 03 EE=EV/TEHP . 

PI=3. 1415926536 
AC13=-1,/SQRT(3»3 
AC23=-AC13 
M0STPS=0 
MSTEP=1000 

SET UP SUMS OVER TRAJECTORIES 

IFCKA.EQ.03 GO TO 250 
J A M A X = 2 
JBM AX=2 
K AM AX=MA 
KBM AX=MB 
NUHBER=MA*HB*4 
IFCNN.E0.13 WRITECH, 9903 


C 

C 

C 


3 O 

■*1 § 
-9 5) 

o 5 
of 

23 £— 

<0 *a 


3, 


m 


990 


938 

C 


HA, HB, NUMBER 

FpRMATC/lX,I4,16H ALPHA-INTERVALS, 3X, 14, 15H BETA-INTER VALS , 6X, 
L 5HHENCE » 14, 35H TRAJECTORIES FOR EACH ENERGY-VALUE 3 

IFCME.EO-03 GO TO 200 

HE2=2*HE 

JEMAX=2 ' 

KEMAX=ME 

IFCNN. 60,13 WR ITEC M, 98 8 3 HE,HE2 
FORMATC1X, 14 ,27H ENERGY INTERVALS 


AND HENCE, 14, 14H ENERGY VALUES3 


GO TO 300 


C 

C SINGLE VALUE OF ENERGY 
C 

200 


996 


JEM AX= 1 
X EM AX= 1 

IFCNN-E0.13 WRITE CH, 9863 EV.EE 
FORMAT C1X, 31H HONOENERGETIC CASE 
DIMENSIONLESS VALUE, E16. 43 


WTIH ENERGY, 1PE16. 4, 30H VOLT<^ CS 




\ 

i 



ISN 

0049 

C 

r 


GO TO 300 




L 

c 

r 

SINGLE TRAJECTORY ONLY 


ISN 

0050 

250 

JAHAX=1 


ISN 

0051 



JBMAX=1 


ISN 

0052 



JE M AX= 1 


TSN 

005 3 



XAHAX=1 


ISN 

0054 



KBMAX=1 


ISN 

0055 



KEMAX=1 


I SN 

0056 



AL=AL1*P 1/180. 


ISN 

0057 



BE=BE1*PI/180- - 


ISN 

0058 



WRITECM* 904 ) AL1 , AL, BE1 , BE , EV , EE 


ISN 

0059 

984 

FORMAT C/1X , 17HSINGLE! TRAJECTORY 






1/1X, 7 IIALPHA =,F20.8 ,12H DEGREES. 0R.F20.8 , 8H 

RADIANS 





2/1X, 7MBET A =,F20.8 ,12H DEGREES, 0R.F20.8 , 8H 

RADIANS 





3/ IX < OMENERGY =, 1PE16.4 , 30H VOLTS. OR 0IHENS10NLESS VALUE , E16.4 3 

ISN 

0060 



SIN A=S IN CAL) 


ISN 

0061 

r 


COS A=COS CAL) 




U 

c 

r 

SUH 

OVER ENERGY, BETA, AND ALPHA 


ISN 

0062 

L 

300 

CONTINUE 


ISN 

0063 



DENST=0. 


ISN 

0064 



DO 1001 KE=1 .KEMAX 


TSN 

0065 



DO 1001 JE=1 , JEMAX 

* 

ISN 

0066 



D ENS=0 - 


ISN 

0067 



NOESC— 0 


ISN 

0068 



DO 1000 KB=1 .KBHAX 


ISN 

0069 



DO 1000 JD=1 , JBMAX 


ISN 

0070 



DO 1000 KA=1 .KAMAX 


ISN 

0071 

r 


DO 1000 JA=1 « JAM AX 




u 

c 

r 

INITIAL POSITION" 


ISN 

0072 

U 


Z=0. 


ISN 

0073 



X=XSA VE 


ISN 

0074 



Y=YSAV£ 


ISN 

0075 

r 


IFCHA.EQ.03 GO. TO 320 


ISN 

Q077 



CA=(AC JA3 +■ FLOAT ( 2*KA - 1 - HA33/FL0ATCHA3 


TSN 

0078 



S IN A=S QRTC - 5 *C I. +CA 3 3 


ISN 

0079 

r 


COS A=S QRTC1- - SINA**23 


ISN 

0080 

u 


CBETA = CA( JB3 + FLOAT C2*KB - 1 - HB33/FLOAT CH83 


ISN 

0081 

c 


BE=PI*C 1 • *■ CBETA3 


ISN 

0082 

c 

320 

XD0T=SINA*C0SCBE3 


ISN 

0083 



YDGT=SINA*SINCBE3 


ISN 

0084 



ZDOT=COSA 


TSN 

0085 



INI =0 


ISN 

0086 

r 


CALL INTERP 


ISN 

0087 



IFCIGOUT.GE.l.AND.IGOUT.LE.IIA.AND-JGOUT.GE-l.AND 
1 AND.XGOUT.GE. 1. AND-.KG0UT.LE-KKA3 GO TO 340 

* JGOUT-LE- JJA . 

ISN 

0089 

330 

WRITE CM, 99993 




1 


; 

UU 'JO 

99*39 

I SN 

0091 


I SN 

0092 


1 SN 

0093 

C 

tsn 

00 94 

340 

ISN 

0095 


ISN 

0096 


ISN 

0097 


ISN 

0098 

C 

ISN 

0100 


ISN 

0101 

C 

ISN 

0102 

350 

ISN 

0103 


ISN 

0104 


ISN 

0106 

C 

ISN 

0107 

400 

ISN 

0109 

C 

ISN 

0110 


ISN 

Oil 1 


I SN 

0112 


ISN 

0113 


I SN 

0114 


ISN 

01 1 5 


ISN 

0116 


ISN 

0117 


tsn 

0118 


I SN 

01 1 9 


ISN 

0120 

C 

C PR 
C 

ISN 

0122 


ISN 

0123 

982 



C 

ISN 

0124 


tsn 

0125 

980 



C 

ISN 

0126 


ISN 

0127 

8 B 8 



C 



C 



C 

I SN 

0128 

490 

ISN 

0130 

500 

ISN 

0131 


ISN 

0132 


1 SN 

0134 


ISN 

0136 


ISN 

0137 

998 


FORMAT C /////IX » 4 3HOME OF THE IG-JG-KG INDICES IS OUT OF RANGE! 
WRITE CM ,000! KST EP , X, Y , Z , XIDOT , YDOT , ZDOT , I GOUT, JG0UI,XGOUT « PHI 
WRIT£CM,982)KE,J£,KB, JB , X A , JA , U E I , AL1 , EV , P VOLTS 
STOP 

INT = 1 

PH ISAV=PHI 
SPEEU=0. 

PH I0LD=PHT 

IFCME.GT.O) GO TO 350 
E = EE 

GO TO 400 

CE=CAC JE) * - FL0ATC2*KE-1-HE))/FL0ATCHE) 

E-Ci-+CE>/tl.**CE) 

IFCXHACH.GT.l.) £=XHACH**2*C1.+CE)/C1.-CE) 

E=£ + AHAXlCPHlt 0.) 

IFCE.LT.PHI) GO TO 1001 
SP£ED=SQRTCErPHl) 

X00T=SPEED*SINA*C0SC8E) 

YD0T=SPEE9*5INA*SIHCBE) 

ZD0T=SPEE0*C0SA 

AL=ARCOSCCOSA) 

AL1=AL*100./PI 
BE1=BE*180./PI 
EV = E* r TEHP 

P VOLTS =P HIS AV^TVOLTS 
ZGL0=Z 
K S T E P = 0 

IFC NPR INT.NE -2.AND -NPRINT. NE-33 GO TO 459 

INT INITIAL CONDITIONS OF TRAJECTORY 

WRITECH, 902) KE , JE , XB, JO ,KA , J A, 8 El , AL1 , EV , PVOLTS 
/-O f }!) A ^<i x *52HXE, JE, K8,JB, KA » JA, BETA, ALPHA, ENERGY, P9TENTIAL= 
i’r H»3CI3,I2),1PE22.8 ,4H DEG,4X,£22.8 ,4H DEG,8X ,E16.4,2H V,4X, 

2 E16-4 S 2H V) 

WRITECH ,980) 

FORMAT ( 9X, 95HSTEPS 

i yoot zoor 


X 

IG 


JG 


PHI) 


XDOT 


F«^ 1 ^f?*?^? 5 !f!fI?i?: 5 :lig^l!T |? T ’ ZDar * rGOOT - JcauT * KGOUT ’ pi ’ 1 

TAXE A STEP 

IF CXSTEP.EO.O) GO TO 550 
CALL ORBIT 
K S T EP = KS TEP + 1 

, IF ,^;ir,? I t , JA, E ,9*^L WRITECM * 808:) XSTEP, X, Y,Z ,XDOT, YOOT, ZOOT, IGOUT, 

l JuiJUi jAuuUi i "nl 

IFCKSTEP.LE.HSTEP) GO TO 550 
WRITECH, 998) HSTEP 

FORMAT C/////1X, 9HHORE THAN, 16,1911 STEPS - HENCE STOP) 


// 


I 


ISN 0138 


STOP 


ISM 0139 
ISN 0141 
ISN 0143 


ISN 

ISN 


ISN 
ISN 
ISN 
I SN 

ISN 

ISN 

TSN 
I SN 
ISN 
ISN 


0145 

0147 


0149 

0150 

0152 

0153 

0155 

0156 

0157 

0159 

0160 
0161 


ISN 0163 


TSN 

ISN 

ISN 

ISN 

ISN 

ISN 


TSN 

ISN 

ISN 

ISN 

ISN 

ISN 


ISN 

ISN 

ISN 

ISN 

ISN 

ISN 


0165 

0167 

0168 
0170 
0172 
0174 


0176 

0170 

0179 

0181 

0103 

0105 


0187 

0109 

0190 

0192 

0194 

0196 


C 

C 

550 


538 


540 

C 


C 

560 


C 

570 


IFCZ-EG.O,-AND-ZDOT.LT.O- 
l-AND-Y-GE-YPClJ-AND. Y.LE.YPCJV33 GO TO 600 

IFCCX.LE-XXCD.AND-ZD0T.LT.0-3.0R. 
ICX-GE-XXCII3-AND.Z00T.LT.0. 33GO TO 600 

IFCCX-LE.XXC13-AND.XDOT.LT.O..ANO.ZOOT.GD OD.OR, 
1 CY. LE.YYC13.AN0.YD0T.lt. 0.3. OR. 

2CX-GE-XXCII3 .AND.XDOT.GE.O.-AND.£OOT.GT.O*3.0R. 

3 CY.GE.YYC JJ3.AND.YD0T.GT.0. 3. OR. 

4 CZ.GE-ZZCKFO.AND.ZDOI.Gl -0.33GD TO 700 


? 


^ 0 


.S' 


ISN 0190 


IFCSKPRFL.EC.13 GO TO 538 

IFCCCY.LE.CY YCH3C3--5*Z3 3.AND.CY„GT.CYYCHBC3-.5* i ZZCXUK3333.0R. 
1CCY.GE.C YYCMBD3+-5*Z3 3-AND.CY.LT„CYYCMBD3 + .5*ZZCXUK3 33 3 
2-AND-X-GE.XXCILX3-AND.X-LE.XXCIUX33 GO TO 600 
CONTINUE 

IF CZ-NE.O..OR .ZDOT.GE.0.3 GO TO 540 
ZDOT=-ZDOr 

IF CNPRINT . EO. 33 HRITECH,8883 KSTEP*X, Y,Z,XDQT,YDflT*ZDOT*IGGUT » 

1 JGOUI »KGDUT »PHI 
GO TO 590 
CONTINUE 

IF CKSTEP-EQ.03 GO TO 500 

PHrOLO=PMI 

CALL INTERP 

IFC IGOUT . LT. 1. OR . I GOUT . GT -IIA.OR-JGOUT . LT . 1. OR. JGOUT. G T. JJA. OR. 
1KG0U T • LT -I.0R.KG0Uf.GT.KKA3 GO TO 330 
IFCNTIHE-LT.1.0R.NTIHE.GT-63 GO TO 580 

IFCNTIHE-NE.1-AND.NTIHE.NE.23 GO TO 560 
XD0TS=XD0T**2 *- PHIOLD-PHI 
IFCXD0TS.EQ.0.3 XDOr=0. 

IFCXD0TS-GT.0.-AND.X00T-NE.0.3 tl>OT=SORTCXDOTS3*SIGNCl.,XDOT3 
IFCXDOTS.LT.O. . AND -XOOT -NE -0 -3 XuOT=-XDOT 

IFCNP RINT.E0.3.AN0-XD0TS.LT. 03 WR I TECH, 8883 KSTEP,X, Y, Z.XDOT, YDOT, 
1 ZDOT.IGUUT, JG0UT,KG0UT,PHI 

IFCNTIHE -NE- 3.AND.NTIHE.NE.43 GO TO 570 
Y00IS=YD0T**2 + PHIOLD-PHI 
IFC YOD IS.E0.0.3 Y001=0- 

IFCYOOTS .GT.O..ANO.YOOT.NE.O-3 YDOT=SQRT CYD0TS3*SIGNC1 . * YD0T3 
TFCYOOTS-LT -O..AND-YOOT.NE.O-3 YOOT=-YDOT 

IFCNPRINT.EO.3-ANO.YOOTS-LT.03 WRITE CH, 808 3 KSTEP ,X, Y, Z, XDOT, YDQT. 
1 ZOUF, IGOUT, JGOUT, KGOUT, PHI 

IFCNTIHE-NE.S.nND.NTlHE.NE.63 GO TO 590 
ZD0TS=ZD0T**2 + PHIOLD-PHI 
IFCZDOTS.E0.0.3 Z001=0. 

IFCZD0TS.5T.0..AN0.ZD0T.NE.0.3 Z00T=S0RTCZ0CTS3*SIGNC1 - ,ZOPT 3 
IFCZDOTS .LT.O..AND.ZOOT.NE.O-3 Z00T=-ZD0T 

IF CNPRINT- EQ.3-AND.ZD0TS.Lr. 03 WRITECH,8383 XSTEP ,X, Y,£ *XDOT, YDOT, 
1 ZDOT, IGOUT, JGOUT, KGOUT, PHI 
GO TO 590 


© o 

**1 -J 



•• I 


I 


TSN 0199 
ISN 0200 
ISN 0201 
ISN 0202 
ISN 0203 

TSN 0204 
ISN 020b 

ISN 0207 

ISN 0209 


ISN 0210 
ISN 0211 
TSN 0213 
ISN 0214 
ISN 0215 


ISN 0216 
ISN 0217 
ISN 0219 
I SN 0221 
ISN 0222 
TSN 0223 
ISN 0224 
•ISN 0225 
ISN 0226 

ISN 0228 
ISN 0229 
ISN 0230 
TSN 0231 
ISN 0233 
TSN 0234 

ISN 0235 
TSN 0237 
ISN 0230 

TSN 0239 
ISN 0240 
ISN 0242 
ISN 0243 
ISN 0244 
ISN 0245 
ISN 0246 
ISN 0247 
ISN 0249 
TSN 0249 


580 WRITECH, 9973 NTIHE 

997 FORMAT C////IX, 17HTRQUBLE - NTIHE =,I3,19H = OUT OF RANGE 1-63 

WRITECH, 8073 KSTEP, X , Y , 2, XQOT, YDQT, ZDOT , IGOUT, JGQUT, KGOUT, PHI, SAVE 
887 FOR MAT C 9X , I 5 , 1P6E1I - 3 , 3 1 6 , El 1 -3 , *S A VE = * , E 1 8 . 103 

STOP 
C 

590 CALL INTERP 

TFCIGGUT.LT.l.OR.IGOUr.GT.lIA.OR.JGOUTiLT.l.QR.JGnUT.GT-JJA.OR. 
1KG0UT.LT.1. OR. KGOUT. GT.KKA) GO TO 330 
IFCNPRINT.EO.33 WRITECH, 8883 KSTEP , X £ Y, Z, XOOT,YDOT, ZDOT# IGOUT, 

1 JGDUT , KGOUT » PHI 
GO TO 500 
C 

C PARTICLE IS A8SQRBED ' 

600 CONTINUE 

IFCNPR INT-NE. 2- AND-NPRINT.NE-33 GO TO 1002 
FATEC13=EN01C13 
FA T 6 C 2 3= ENO 1 C 23 
GO TO 750 
C 

C PARTICLE ESCAPES 
C 

700 CONTINUE 

IFCNPRINT-E0.13 GO TO 720 

IFCNPRINT.NE-2-AND.9PRINT.NE. 33 GO TO 740 
FATECJ 3=EMQ2C13 
FATEC23=E9D2C23 
GO TO 74 0 

720 UR I TtCK, 9B23 KE , JE . KB. JB , KA, JA , BEl» AL1 , EV, PVOLTS 
740 NOESONOESC *• 1 

IFCM6-EQ.Q3 GO TO 750 
C 

CSANGL=ZnBT/SORTCXDOT**2+YDOT**2+ZDOT**23 
XPDN=-2.*XHACH*SQRTCE3*CSANGL - E - XHACH**2 
CO E FA = SP F;E 0**2 /FLOAT C NUMBER 3 
IFCA B SCX P0N3-GT. 36-3 GO TO 1000 
ADO =C0EFA*EXPCXP0N3 
DENS=DENS + ADO 
C 

750 IFCNPRlNT.Nt-2-AN0.NPRINT.NE. 33 GO TO 1002 

WRITE CM, 089 3 FATE, KSTEF, X, Y, Z, XOOT ,YDOT, ZDOT, IGOUf , JGOUT, XGOUT, PHI 
889 FORMAT C IX, 2A4, 15, 1P6E1I. 3 , 316 , E 11 . 33 
C 

1002 CONTINUE 

IFCM0STPS.GE.KSTEP3 GO TO 1000 

KES-XE 

JES=JE 

K B S =K 8 

JO S = JB 

KAS=KA 

J AS =JA 

MOSTPS=KSTEP 
1000 CONTINUE 
C 

C END OF SUH OVER ANGLES 
C 


ISN 0250 


FRACT=FLOATCNOESC3/FLOATCNUHBER3 


TSN 

0251 


ISN 

0252 

970 



£ 

TSN 

0253 


1 SN 

0255 


ISN 

0257 

976 



C 

TSN 

0250 

800 

TSN 

02 6 0 


ISN 

026 L 


ISN 

0263 


ISN 

0266 

1001 

C 

ISN 

0265 

C 



C 



C 

TSN 

0267 


TSN 

0268 

972 

ISN 

0269 


L SN 

0270 

976 

TSN 

0271 


IS N 

0272 



WRITECH, 978) NOESC , NUMBER ,FRACT, EV.DENS 

FnRHArC/lX,16!IRATI0 ESCAPING =,.15* 71! DU. OS'rl5.16H OR A FRACTION, 
1 Fl3.fi, 1411 AT ENERGY E =, F13.il, 6H VOLTS* X * 6HCQ£NS=, I PE16 . 6, 1IO ) 

IFCHPRINT-EO-O) GO TO flOO 
IF C HE. N 6.0) WRITECH, 976) 

FORMAT C IXi BOUDENS IS THE SUM OF AOD=SP EE0**2*EXP CXP0N)/fJUH3ER OVER 
1 A HEMISPHERE OF DIRECTIONS//) 


IFCHE.EQ.Q) GO TO 1001 
C0£FE=2./C1. - CE)** 2/ FLOAT CHE) 

IF C XM ACH . GT . 1 . ) CQ£FE=CQEFE*XHACII**2 

DEN ST =OENST + C0EFE*OENS 

CONTINUE 

IFCME.EO.O) DENST=SPEED**2*FRACT 

TRAJECTORY WITH HOST STEPS- PRINT X AND J INDICES. 


WRITECH, 972) HOSTPS . XES , JES ,KBS , JOS , XAS, JAS 
FORHATC///IX, T5,3CI3,I2),29H =HOSTPS, KE,JE, KB.JB, KA , JA) 

WR I TEC H, 970 XS AVE , YSAVE , PllISAV , DENS T , PAR TCL 

FURMATC/IX, 2 6HAT DIMENSIONLESS XtY.PHI = , 3F1 2.6, 111, , 5X . 1PE16. 
1 3 311 = NORMALIZED CURRENT DENSITY F0R*2A5//) 

RETURN 

END 


If « 


l 


V uuil It J 


05/360 FORTRAN H 


DATE 


C 



ISN 

c 

TSN 

1 

o .. 



ISN 

ISN 

1 . 

o 


o 

ISN 

ISN 

ISN 

ISN 

o 

ISN 

ISN 

ISN 

c 

ISN 

ISN 

ISN 

< , 
(! 

ISN 

ISN 

ISN 

ISN 

« -1 

ISN 

ISN 

ISN 

ISN 

TSN 

ISN 

ISN 

* * 

ISN 

ISN 

O .* ! 


0 


o 

ISN 

ISN 

ISN 

1 

ISN 

ISN 

ISN 

ISN 


COHPILER 


0002 



C 


C 


C 

0003 


0004 


0005 



C 

0006 


000 7 


0008 


0009 



c 

0010 


0 0 L 1 


0012 



c 


c 


c 

0013 


0015 


001 7 



c 

00] 9 


0020 


0021 


0023 

10 


c 

0024 

100 

0026 


0028 

101 

0029 


0031 


0032 

102 

0033 



c 

0035 

103 

0036 

104 


c 


c 


C 


C 


C 


c 

0037 


0039 


0041 



c 

0043 


00 4 4 


004 5 


004 7 

20 


C 


OPTIONS - NAME = MAIN, 0PT= 02, LINECNT=60 , SI 2E=0000K, 

SOURCE, EBCDIC, NOLIST, NO OECK, LOAD, HAP, NOEDIT, NOID.NOXREF 
SUBROUTINE INTERP 

INTERPOLATION WITHIN GRID 

COMMON7BK/IIMpIIP,JJH s JJP,KK»NTOT , IV,JV,IT »JJ,M»N»VPC303» 
1XYZC208G.3D, WC 3 0 , 20 , 1 0 ) , XP C3 03 , XMC 10) , Y P C203 , YMC10:> , ZZC 1 03 , 
2XXC403 , YYC303 ,ILX, IUX , KUK , MBC , MB D , VRF , NFPS , SKPRFL , SKP LST 
COMMON /ORB/XDOT, YDO I , ZDOT , XI , X2 • YI , Y2 , Z1 , Z2 , X , Y , Z , PHI , NTIHE ,S A VE 
COMMON /INTER /INI ,IIA , JJA ,KKA,IG0UT, JGOUT , K GOUT, XA , YA , ZA , 

1XIC30) «YJC203»ZKC103 

IG0UT=0 
JGOUT = 0 
KGOUT= 0 
NCH=0 

X A=X 
YA = Y 
ZA=Z 


LOCATE XA 

IFCXA.EQ.XICIIA33 
IFCXA. EQ.XICIIA33 
IFC IM7.NE.03 

DO 10 1=2, IIA 
I G=I-1 

IFCXA.LT. XICI33 
CONTINUE 


IG=IIA-1 
GO TO 103 
GO TO 100 


GO TO 103 


IFCXA. GE.XICIG+133 GO TO 102 ’ 3 

IFCXA. GE.XICIG33 GO TO 104 
IG=IG- 1 

IFCXA.LT. XICIG33 GO TO 101 

GO TO 103 

IG=IG+1 

I F CXA.GE.XICIG<-133 GO TO 102 

NCH = 1 
CONTINUE 

ACCEPT IF XICIG3 LESS THAN DR EQUAL TO XA LESS THAN XIClG+13. 
LOCATE YA 

IFCYA. EQ.YJCJJA)) JG=JJA-1 
IFCYA.EQ.YJCJJA33 GO TO 203 
IFCINT.NE.03 GO TO 200 


O o 

^ 33 
O $ 
i— 

S -o 

*8 

3® 


DO 20 J=2, JJA 
JG= J-l 

IFCYA.LT . YJC J3 3 
CONTINUE 


GO TO 203 


'• i 


/3 


I 



« 


ISM 

005 0 

200 

IFCTA-GE.YJCJG+133 GO TO 202 







I S N 

0050 


IFCYA.GE .YJC JG3) GO TO 205 







I SN 

0052 

201 

JG = JG - 1 







1 SN 

00 53 


1FCYA.LT. YJC JGJ3 GQ TU 201 







I SN 

0055 


GO TO 203 







1 S M 

005 6 

202 

JG = JG *■ i 







ISM 

0057 

p 

IFCYA.GE. YJC JG+133 GO TO 202 







l S N 

0059 

203 

NC 11 = 1 







I SN 

005 0 

205 

r 

CONTINUE 









u 

c 

c 

ACCEPT IF YJCJG3 LESS THAN OR 

EQUAL 

TO 

YA 

LESS 

THAN 

YJCJG+13 



c 

c 

p 

LOCATE ZA 







ISN 

0061 


IFCZA.EO-ZXCXXA);) KG-KK A-l 







ISM 

0063 


IFCZA. EO.ZKCKKA33 GO TO 303 







ISN 

0065 

p 

IFCINT.NE.G) GO TO 300 







ISN 

0067 


00 30 K=2,XXA 







I SN 

0068 


KG=K-1 







I SN 

0069 


IFCZA. LT.ZXCK33 GO TO 303 







I SN 

0071 

30 

p 

CONTINUE 







ISN 

0072 

300 

IFCZA.GE.ZKCKG+133 GO TD 302 







I SN 

00 7 5 


IFCZA. GE .ZKCKG33 GQ TO 305 







I SN 

007 6 

301 

KG=KG-1 







ISN 

0077 


IFCZA.LT. ZKCKG33 GO TO 301 







ISN 

0079 


GO TO 303 







TSN 

0080 

302 

KG = KGH 







ISN 

0081 

r 

IFCZA.GE .ZXCKG+13) GO TO 302 







ISN 

0003 

303 

H C 1 1 = 1 







I SN 

00 85 

305 

r 

CONTINUE 









L> 

c 

r 

ACCEPT IF ZXCKG3 LESS THAN OR 

EQUAL 

• 

TO 

ZA 

LESS 

THAN 

ZKCKG*-!) 



u 

c 

p 

LOCATE LIVE AND BOX 







TSN 

0085 

Lr 

X1=XI C 1G ) 







I SN 

0006 


Y1=YJC JG 1 







I SN 

0087 


Z1=ZKCKG3 







TSN 

0080 


x2=xicigh n 







rSN 

0089 


Y2=YJC JG + n 







ISN 

0090 

P 

Z2 = ZKCKG+n 







ISM 

0091 

U 

IFCX.NE.Xl.OR.XDOT.GE.O.3 GO 

TO 400 






ISN 

0093 


IG=IG-1 







I SN 

0095 


X 2 = X 1 







I S N 

0095 

p 

X1=XI C IGJ 




- 



1 S N 

0096 

500 

IFCY.NE. YI.OR.YDOT.GE. 0. 3 GO 

TO 500 






ISN 

0098 


JG=JG-1 







ISN 

0099 


Y2= Y1 







ISN 

0100 

p 

Yl= YJC JG3 







ISN 

0101 
III 1 i 

500 

IFCZ.NE.Z1.0R.ZD0T.GE.0.3 GO 

TO 600' 







/< 


1 3N 

tsn 
r sn 

tsn 

I SN 
ISN 
TSN 
ISN 
ISN 


UIUJ 


KG=KG-1 

0104 


Z2=Z1 

0 L 0 5 

r 

Z1=ZK(KG) 

0106 

600 

PHI=VVCIG, JG.KGO 

0107 


IGOUT = IG 

0100 


' JGOUT=JG 

0109 


KGOUT=KG 

0110 


RETURN 

0111 


END 



DA 


COMPILER OPTIONS - NAHE = M A I N , 0 PT = 02 , LINEC NT=60 , S I Z E=G 000K * 

SOUR CE.EBCDIC.NOL 1ST, NOD ECK, LOAD, HAP, NOEDIT, NO ID, NOXREF 
SUBROUTINE POWER 

CURRENT DENSITIES AND POWER LDSS 

COMMON/BK/IIK,IIP, JJH, JJP ,KK , NTOT, IV, J V , II , J J,H , N, VP C3 0) , 
lXYZC2080,3)fVVC30i20,10),XPC30),XHC10),YPC20)tYHClQ),ZZClO)» 
2XXC40) ,YYC30),ILX«IUX,KUK,MBC,MBD,VRF»NFPS, SKPRFL , SKPLST 
COM MON/C P/NP ft IHT.NPTS.MA , MB , ME , KMAX , X P T , YP T, AL1 , BE1 , E V , SHACH , 

1 TVOLTS.CUR, XMETER 

COMMON/CD/P VOLTS .XHACM .DENST, NN. PARTCL C2) , PART1C2) .PART2C2) 
DIMENSION A C 2) 

IFCNPTS.EQ.O.OR.MA.EQ.O) WRITECH , 997 ) X P I' , YPT , ALI , BEl , EV 
FORMAT (/ IX , 9HX AND Y = , 2 Fl 0„5 , 20X , 19 HA LP HA , BETA , ENERGY =,3F20.5) 
IFCNPR INT.EQ.O) WRITECH, 990) 

IFCNPRINT.EO.l) WRITECH, 991) . ... 

IFCNPRINT.E0.2) WRITECM, 992) 

1FCNPRINT.EQ.3) WRITECM, 993) 

FORMAT C/3BH NPRINT=0 MEANS NO 
FORMAT C/53 H NPRINT=1 PRINT 
FORMAT C/56H NPRINT=2 PRINT 
IS) 

FORMAT C/52H NPRINT=3 MEANS PRINT EVERY STEP OF ALL TRAJECTORIES) 
IFCTVOLTS.EQ.O.) RETURN 
XMACH=SH ACM . 

NON-DIMENSIONALIZE THE POTENTIAL DISTRIBUTION. THEN RESTORE AT END. 

DO 200 K=1,KK 
DO 200 J=1,JJ 
DO 200 1=1,11 

VVCI.J.K )=VVCI, J,K)/TVOLTS 
CONTINUE 
C 

C DEFINE THE PANEL POINTS AT WHICH THE CURRENT AND POWER IS EVALUATED 

C CASE OF A SINGLE POINT 
C 


PANEL SUB^AREAS 


ISN 

0002 

C 

c c 
c 

isn 

0003 


I SN 

0004 


TSN 

0005 


isn 

0006 


ISN 

0007 


ISN 

0009 

997 

ISN 

0010 


ISN 

0012 


ISN 

0014 


ISN 

0016 


ISN 

0018 

990 

ISN 

0019 

991 

ISN 

0020 

992 

ISN 

0021 

9 93 
C 

ISN 

0022 

C 

ISN 

0024 

c 

C N 
C 

ISN 

0025 


ISN 

0026 


ISN 

0027 


ISN 

0028 


ISN 

0029- 

200 


TRAJECTORY PRINTING) 

INDICES OF ESCAPING TRAJECTORIES ONLY) 
FIRST AND LAST STEPS OF ALL TRA JECTORIE 


ISN 

0030 

r 

IFCNPTS.EQ.0.OR.MA.EQ.O) 

COEFM = XMETER**; 



u 

c 

r 

CASE OF MULTIPLE. POINTS FOR 

INTEGRATION OVER 1 

ISN 

0032 

u 

JVH = 1 


ISN 

0033 


I VM = 1 

* 

ISN 

0034 


IFCJV.GT.l) JVN=JV-1 

O 

I SN 

0036 


IFC1V.GT.1) IVM=IV-1 


ISN 

0033 


NA = 0 


ISN 

0039 


NAREAS=IVM*JVM 

O 3? 

TSN 

004 0 


I POWER=0 - 

pf 

ISN 

004 1 


T C U RNT = Q . 


ISN 

0042 


T AR E A = 0. 


ISN 

00 43 


NN=0 


ISN 

004 4 


00 500 J = 1 , J V H 

£ © 

ISN 

004 5 


DO 500 1=1, IVH 


ISN 

004 6 


NA=NA+1 

^19 






.i£tj._ 


*ZiilZd.j£ 


rsN 

00 4 7 


ISN 

0048 


ISN 

005 0 


ISN 

0051 


ISN 

0052 


ISN 

0053 


ISN 

005 4 


ISN 

0055 

250 

ISN 

0056 


ISN 

0057 


ISN 

00 5 8 


ISN 

Q Q 5 9 
006O 


ISN 

26G 

ISN 

006 1 


ISN 

0062 

270 

ISN 

0063 


ISN 

00 6 4 


TSN 

00 6 5 


ISN 

0066 


ISN 

0067 


I SN 

0068 


ISN 

0069 


rsN 

007 1 


ISN 

0072 


ISN 

0073 


ISN 

0074 


TSN 

0075 


ISN 

0076 


ISN 

0077 

C 

C 

c c 
c E' 

ISN 

0078 

300 

ISN 

0079 


ISN 

0080 


ISN 

008 1 


ISN 

0083 


TSN 

0084 


ISN 

0085 


ISN 

0086 


ISN 

0087 


ISN 

0088 

C 

ISN 

0089 

995 



C 

ISN 

0090 


ISN 

0092 

99 4 



C 

ISN 

0093 


TSN 

0094 


ISN 

0095 

988 

ISN 

0096 



NP = 0 

IFCNPTS.EQ.O.OR.HA.EQ.Q9 GO TO 250 
PO = 0 
CU = 0 

AC13=-l./SQRTC3-3 

A(23=-AC1) 

GO TO 260 
CONTINUE 
JXM AX= 1 
J Y H A X = 1 
KHAX=1 
GO TO 270 
J X H A X = 2 
J Y M A X = 2 
CONTINUE 

DO 400 K Y=l» KHAX . 

DO 400 KX=1 , K MAX 
D'O 4 0 0 JY = L * J YHAX 
DO 400 JX=1 , JXHAX 
NP = NPH 
N N = N N + 1 

IFCNPTS.EO.0.QR.HA.EQ.O3 GO TO 300 

CX=CAC JX3 + FLOAT C2 fc KX - I - KHAX 03 /FLOAT CKHAX9 

CY = CACJY3 + FLOAT C 2*K Y - 1 - KM AX33/FL0 AT C K MAX3 

XPT = CXPCHT3«XPCI)3/2.*CX «- CXP C1+13+XPCI33/2. 

YPT = (YPC J+13-YPCJ33/2.*CY * C YP C J U3 *YP C J33/2. 

COE F = CXPCHT3-XPCI33*CYPCJ+13-YPCJ33 

AREA = COEF * XHETER**2 

COEFH = ARE A /4 - /FLOAT C KHAX **23 


Dtxr 


*+■ V 


LT7ZJ~ ^ Nil /(/lA*n 


•UTE EACH CURRENT DENSITY AND KU'LTIPLY BY LOCAL POTENTIAL TO 
-UATE POWER DENSITY 
CALL DEN 

DENCUR=OENST*CUR 
POWDEN=P VOLTS *DENCUR 
IF C MA« EQ .03 GO TO 600 
XPTM=XPT*XHETER 
YPTH=YPT*'X METER 
XPM=XPCI3*XMETER 
XPPM=XPCIH3*XMETER 
YPH=YPCJ3*XMETER 
YPPM=YPC J<T3*XHETER 




■t 

3^ 


— NNy/'lJlfA. 


FORMAT C6X,12HAT POINT NQ.,I3,10H, WITH X . = , F10. 5 , 130 HETERS, 
1.F10. 5*2711 METERS, ANO COEFFICIENT =,F10.5,14H SQUARE ME TERS3 

IFCNPTS.GT.0.AND.MA.GT.03 URITECH.9943 NA,XPM 
FORMAT C /5X.16H IN SUB-AREA "" 

1 FIT. 5, 30 T0,F10.5, 9H3 METERS ,, 4X , 13HAN0 BY Y 

2 FIT). 5, 30 T0.F10.5, 8H3 METERS3 


Y = 


jLickni^s-u MMjArn,XPPH,YPMfYPPH 
HO. ,13, IX.17HB0UNDED ‘ 


BY X IN C, 
IN C, 


,XPTM, YPTH, COEFH 

)WDEN,PARTCL 


WRITE CM ,9953 NP 
WRI-TECM, 9883 P VOLTS , DENCUR , POl 

FORMAT C6X.53HTHE VOLTAGE, CURRENT DENSITY, AND POWER 
1/6X»1PE16. 4, 6H VOL TS,4X,E16. 4,230 AHP/CSQ-HETER3 , 

2 240 WATT/CSQ-METER3, F0R,2A5//3 
IFCNPTS.EQ.03 GO TO 600 


DENSITY ARE = 
AND , E16.4 , 


/f 






If 




4 


4 




J,* 


it 


I SN 

009U 


ISN 

009 9 


ISN 

0100 

400 

ISN 

0101 


TSN 

0102 

C 

I SN 

01 03 


ISN 

0104 


ISN 

0105 

984 

ISN 

0106 

986 



C 

ISN 

0107 


ISN 

0108 


r sn 

0109 


1 SN 

0110 

500 



C 

ISN 

0111 


ISN 

0112 

982 

TSN 

0113 


TSN 

01 14 


ISN 

0115 


ISN 

0116 

980 


c 

c 

c 

c 


ISN 

0117 

600 

ISN 

0118 


TSN 

0119 


ISN 

0120 


TSN 

012 1 


ISN 

0122 

700 

ISN 

01 23 


ISN 

0124 



CU = CU + COEFM*DENCUR 

PO = PO + COEFM * POWDEN 

CONTINUE 

AVCD = CU/AREA 

A V P 0 = PO/AREA I 


WRITECM,986) NA,CU,P 
WRITECM,9843 NA , AREA 
FORMAT C IX, 18 H IN SU 

1 METERS, /52H THE AVE 

2 E16.4,19H AHP/CSQ-M 
FOR MAT C /1X.18HIN SU 

1 =,1PE16.4,12H AMP, 


O.PARTCL 
,AVCD,AVPD 
Q-AREA NUMBER, 13, 
RACE CURRENT DENS 
ETER3 AND,E16.4,1 
B-AREA NUMBER, 13, 
AND,E16.4 ,14H 


8 H OF AREA, 1PE16-4, 15H SOUARE 
ITY AND POWER DENSITY ARE =, 
6H WATT/CSQ-METER33 
28M THE CURRENT AND POWER ARE 
WATTS, FOR , 2A5 3 


TAREA=TAREA + AREA 
TCURNT = TCURNT J- CU 
TPOWER = TPOWER + PO 
CONTINUE 


POWER LOSS 
FDR , 2A 53 


ARE =, 1PE16. 4, 


WRITE CM, 9823 TCURNT , TPOWER , PART CL 
F0RMATC///1X,34HT0TAL CURRENT AND 

1 12H AMP, AND,E16.4,13H WATT, 

AVCD =T CUR MT/TAREA 
AVPD=TPOWER/TAREA 

WRITECM, 9803 TAREA , A VCD , AVPD 

FORMAT C/1X»26HWITH A TOTAL PANEL AREA OF, 1 PE16. 4 , 15H SOUARE METERS 
1, /IX, 51HTHE AVERAGE CURRENT DENSITY AND POWER DENSITY ARE =, 

2 El 6 • 4 .» 1 9H AMP/CS0-METER3 AND , El 6- 4 , 1 6H WATT/ CSQ-HETER 33 


RESTORE POTENTIAL DISTRIBUTION 


TO DIHENSIQNAL VALUES 


CONTINUE 
DO 700 K=1,KK 
DO 700 J=1,JJ 
DO 700 1=1,11 

VVCI , J,K3=VVCI, J,K3*T VOLTS 

CONTINUE 

RETURN 

END 


IP 


4 


Z. o 


: » 


COMPILER OPTIONS 


ISN 0002 
ISN 0003 


ISN 0004 
ISN 0005 
TSN 0006 
ISN 0007 
ISN 0008 
ISN 0009 
TSN 0011 
ISN 0013 
ISN 0014 
ISN 0015 
ISN 0017 
ISN 0019 
ISN 0021 
TSN 0023 
ISN 0024 
TSN 0025 
ISN 0026 
ISN 0027 
TSN 0028 
ISN 0029 
TSN 0030 
TSN 0031 
ISN 0032 


200 

300 

400 

1000 

450 

3000 

500 


- NAHE= MAIN, 0PT=02,LINECNT=60, SIZE=OOOOK, 

rMOnn„-r T r:9 U ^P?^?99 I 9xC JQLISr * N00ECK » L0AD » HAP » N0 E DIT * N0ID » N 0XREF 
SUURUUTINh LI5TCLST,IP3 

COMMON/OK/IIM, IIP , JJM, JJP ,KK,NTQT,IV, JV.TI, J J,H, N, VPC303 , 
1XYZC20 00,33,VVC30,2G,1Q3,XPC3 03,> XHC 1 03 , Y P C 203 , YHC 1 01 . ZZC1 03 , 
2XXC40) , YYC303 ,ILX, I UX , KUK , MO C , Mil D , VRF , NFP S , SKPRFL.SKPLST 
COMMON/FLD/XC2080.2) , C 0 E FC 2 0 80 , 73 , IND X C 20 8 0 , 63 , SKPCQ 
DIMENSION K0UIC5 3 ,X0UTC5 3 , YOUTC5 3 »ZQUTC53 
DO 500 LINE=1 , 60 
DO 200 NP=1 , 5 

KP=LINE + C.VP-n*60 + CIP-13 * 300 
IFCKP. GT.NTOT. AND- NP. EQ-13 RETURN 
IFCKP «GT . NT0T3 GO TO 300 
NHAX=NP 
K0UTCNP3 = XP 

IFCLST. EQ. 13 X0UTCNP3 = XCKP, 13 
IFCLST -EQ. 23 X0UTCNP3 =XYZCKP,13 
IFCLST - ED - 23 YOUTC NP 3=X YZ CKP , 23 
IFCLST - EQ- 23 Z0UTCNP3 =XYZCKP,33 
CONTINUE 

GO TO C400,4503,LST 

NP = l - NH4 » 

GO TO 500 

^n^^^S°^^ X 9yi CNP) ’ XD1JTCNP:)jYDljT(:NP:> » z OUTCNP3,NP=l,NHAX3 

lUKnAI C 3 C I a i 3r 6. 2)3 

CONTINUE 

RETURN 

END 


i uni RBB 11 


OAT 




ISN 

ISN 


ISN 
ISN 
TSN 
ISN 
I SM 
1 SN 
TSN 
ISN 
I SN 
ISN 
ISN 
ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 
ISN 
1 SN 
ISN 
TSN 
ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 

ISN 


0002 

0003 


COMPILER 

C 


0004 

0005 

0006 

0007 

0008 

0009 

0010 
0011 
0012 

0013 

0014 
0016 

0017 

ooia 

0019 

0020 
0021 
0022 

0023 

0024 

0025 

0026 
0027 
002 0 

0029 

0030 
0032 
0034 
0036 
0038 
0040 

0042 

0043 

0044 
0046 

0048 

0049 
0051 

0053 

0054 

0055 

0057 

0058 


C, 


OPTIONS - NAHE= HAIN, 0PT=02, LINECNT=60, SIZE=0000K, 

SQURCE,EBCDIC,NDLIST,NODECK,LOAD,MAP,NOEOIT,NOID,NOXREF 
SUBROUTINE RELAX 

POINT-SUCCESSIVE OVERRELAXATION METHOD 

COMMON /BK/IIM, IIP, JJM, J J P , KK , .NTOT , I V , JV,II, J J, M, N , VP (3 03 , 

1XY Z ( 20 80 » 33 » VVE 30, 20 , 103 , XP (303 * XH( 103 « YPC 203 , YH(1Q 3 » ZZC103* 

2XX ( 403 , Y Y ( 303 * I LX , IUX * KUK » M8C * HBD, VRF* NFP S • SKPRFL, 5KPLST 
LUMMON/FLD/X (208 0, 23 , C0EFC2080* 73 * 1NDXC2080, 63 , SKPCO 
n m i- r. a = i q 


200 


ISN 0059 


500 

888 


C 

600 

C 

C 


0MEGA=1.9 
EPS = l.E-3 
ITHAX=2000 
I TR = 0 
I PR OLD=0 
I G 0 - 1 
ITR=ITR+1 
DELTA M=0. 

DO 500 N=1 , NTOT 

TF(X(N .23.EQ.13GO TO 500 

X1=X(N,13 

FN=COEFCN, 13 /COE FEN, 73 
FS=C0EF(N,23/C0EF(N,73 
FE=CG£F(N,3)/C0EF(N,73 
FN = C0EF(N,43/C0EF(N, 73 
FU=C0EF(N,.53/C0EF(N, 7 3 
FD=COEF(N,63/COEF(N,73 

NN=INDX(N, 13 
NS = IND X ( N, 23 
N E = I N 0 X ( N . 3 3 
NW=INOX(N, 43 
NU=INDX(N,53 
ND= IND X( N, 63 


© © 
© 

o > 

33 f” 

•O ”3 
CL 

3s 53 
r* m 

3® 


SUM=0. 

IF(NN.GT.03 SUM 
IFCNS-GT-03 SUM 
IF(NE. GT.Q3 SUM 
IFCNW.GT .03 SUH 
IF(ND.GT.O) SUM 
IF(NU.GT.03 SUM 


SUM+FN*X(NN,13 
S'JM + FS*X(..IS,13 
SUM+FE*X(NE,1) 
SUK*-FW*X(NW,1) 
SU M+FD *X( ND ,13 
SUM+-FU*X(NU,13 


X (N , 13 = OMEGA*SUM+(l.-0MEGA3*Xl 
DELTA = ABS(X(N, 13-X13 

IF(ABS(X13-GT.1.E-103 DELTA=ABSC(XCN, 13-X13/X13 

IF(0ELTA .GT. DELTAH3 OELTAH=DELTA 

CONTINUE 

IFCITR.GT.ITMAX) HRITEC M , 8 8 88 3 ITR 
IFEITR.GT.ITMAX3 GU TO 700 
8 FORMAT E////I0H MORE THAN, I4,11H ITERATI0NS3 
IP R-ITR/500 

IFCIPR.LE.IPROLD3 GO TO 600 
IPROLD=I PR 
GO TO 800 

IFCDELTAM.GT -EPS 3 GO TO 200 
ITERATION FINISHED. PRINT AND EXIT. 




r* i* l 


I 




ISN 0061 700 I G (3 =2. 

ISN 0062 000 NFPP=C NT0T/300} +1 

ISN 3063 00 900 IP=1,NFPP 

rSN Cu6 4 WRITECH, 77773 I TR , EPS , D ELT AH , OHEG A 

IS N 0065 7777 FURHA T C i 5H1S GLUT ION AFTER, 16, 2X, 25HITERATI0NS WITH TOLERANCE* 

I Fl2„0,0X,l8HHAXIHUM 01 FFERENCE , FI 2 . 8 , 8X » 6H0MEGA= , F8 . 53 
ISN 0066 CALL LISTCItIP) 

ISN 0067 900 CONTINUE 

J* 

ISN 0068 -&0-T 0 -E 6 0 0,-1- 0 0 0)tT G 0 [F C l (» 0 . £ < 1 ) Or f 0 

ISN 0069 1000 RETURN ^ 

ISH 0070 END 


d 


t 


y+ 


2 * 


v. »i 
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0S2360 FORTRAN H 


DATE 


•«r* 


TSN 

ISN 


ISN 
I SN 
ISM 
S SN 
I SN 
TSN 
ISN 
ISN 
I SN 


0002 

0003 


CQKPILER OPTIONS - NAME= HAIN*0PT = 02» LIMECNT=60 , 5IZE=QGG0iC, 

SOURCE,EBCOIC,NOLIST,NODEC!C, LOAD, HAP , NOEDIT, NO ID, NOXREF 
SUOROUTINE FINDCI » J» JO 


0004 

0005 

0006 
0003 

oooo- 

0010 

0012 

0013 

0014 


COMMON/OK/IIM, IIP, JJH, JJP,KK,MTOT, IV, JV , II . J J, H, N , VP 1 30) , 
C20BO*3),VVC30»2Qil03»XPC3GI)»XH(lu3,YPC203,YHClO}»Z.Z(lO)* 
C 40) i YYC3Q) , 1LX , IUX , KUK , HBC , HBO , VRF , NFPS , SKPRFL ,SXPLS T 
. 1 . 1 = T T * . 1.1 


IX YZ 

2XX. 

IT J J=I I* JJ 
K=N /II JJ* 1 

IF C K .GE. 2 . AND . MDDCN.IIJJ) 
NK I J = N - IIJJ*CX-n 
J = NKIJ/IIH 

IFCJ .GE. 2 -AND. HOOCNXIJ.II} 

I = N K I J - II * C J-13 

RETURN 

END 


. EQ . 0) K=K-1 


- EQ. C3 J=J-1 




i 


»i ii i 


2 V 


*. v. « i. 5* 


u.u JllH ii } 


QS/360 FORTRAN H 


ISN 

ISN 


TSN 

ISN 


ISN 

ISN 

ISN 

ISN 

ISN 

ISN 


ISN 
ISN 
ISN 
ISN 
ISN 
ISN 
ISN 
ISN 
ISN 
I SN 
ISN 
ISN 


0002 

0003 


0004 

0005 


0006 

0007 

0000 

0009 

0010 
0011 
0012 


COMPILER OPTIONS - NAHE = HA IN, 0 P T=02 , LINECNT= 60 , S IZE-OOOOK . 

SOUR CE, EBCD IC, NOL 1ST, NOD ECX, LOAD, HAP, NOEDIT, NOID, NO XREF 
SUBROUTINE ARRAY 

COMMON /OK / 1 1 M* IIP , JJH, JJP ,KK. ,NTOT , I V , J V * II , JJ,H»N»VPC30) , 

1XY a 20 00 ,33, VVC30.20 ,10) ,XPC30) ,XHC10) ,YPC20),YHC10) ,ZZ(10) , 

2XX ( 4 0) , Y Y(30 ) , ILX « IUX » KUK. » HBiC ,MBD ,VRF»NFPS ,SKPRFL»5KPLST 
U)MH0N/FLQ/XC2030,2) , COE FC 2080, 7 ) , IN0XC2080 , 6), SKPCO 
COMMON /CCC/CM , CS , CE » CU, CU , CO, CC , NN, NS , NE « NW « NU,NO 


C 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 


COEFFICIENT ARRAY = COEFCN, 7), WHERE 
COEFCN, 1)=CN CNDRTH=«-Y NEIGHBOR) 
CSOUTH=-Y NEIGHBOR) 

C " E AS T =+X NEIGHBOR) 

C WEST=-X NEIGHBOR) 

C UP=+Z NEIGHBOR) 

C DOWN=-Z NEIGHBOR) 

C = CENTRAL POINT) 


COEFCN, 2) =CS 
COEFCN , 3)=CE 
COEFCN, 4)=CW 
COEFCN, 5)=CU 
COEFCN, 6)=CD 
COEFCN, 7)=CC 


SAVE COEFFICIENTS AND INDICES 


COEFCN ,1)=CN 
COEFCN, 2)=CS 
COEFCN, 3)=CE 
COEFCN ,4 )=CW 
COEFCN, 5)=CU 
COEFCN, 6)=CD 
COEFCN, 7)=CC 


001 3 

0014 

0015 

0016 

0017 

0018 
0019* 
0021 
0022 

0023 

0024 

0025 


INDX C N , 1 )=NN 
INDXCN ,2 )=NS 
I NOXC N , 3 ) = NE 
INDXCN, 4 )=NW 
INDXCN, 5)=NU 
INDXCN, 6)=ND 
IFCSKPCO.EQ- 1) GO TO 20 

, nnrt ND.CD.NS.CS.NU.CU, N , CC, NE , CE , NN ,CN,NU, CU 

1000 FORMAT C/7C1X , Uf C , 14, 2H) = , 1PE10 . 4 ) ) 

20 CONTINUE 
RETURN 
END 


1 


US/360 FORTRAN It 


UATC 


V* 






U* 


*1 0 


<1? 


& , 


COMPILER OPTIGNS - NAME = HAI N, OP T = 02 . LINECNT=60 , SI ZE=00Q0K, 

tsn 0002 su nR o U nrc5i?^^ > t « HOLIST ’ NDOECK,LOM,HAP>HOEOn ' fiOIO,NOXREF 

ISN 000 3 COMMON /OK/ 1 1 H, IIP , JJH, J JP , KK ,NTOT , IV, JV, II , J J,H, N, VP C30) , 

1XYZC2000.3) , WC30,20,10),XPC30) t XHC10) ,YPC20),YHC10) ,ZZCIQ), 
2XXC40), YYC30), ILX, IUX, KUK , HOC, HOD »VRF, NFPS , SKPRFL , SKPLST 
TSN 00 0 A CUMMON/CCC/CN.CS ,CE,CW,CU,CO,CCtNN,NS , NE ,NW,NU,ND 

TSN 0005 NFCIX,JX,KX) = IX+II*CJX-1) + II*JJ*CKX-1) 

ISN 0006 A = 0. 

ISN 0007 C = 0. 

ISN 0000 CALL FINOCI,J,K) 

ISN 0009 II- C I -EO. 13 GO TO 100 

TSN 0011 IFCT -EO. II) GO TO 200 

ISN 0013 NU=NFCI+1, J,K) * 

ISN 0016 NL=NFCI-1,J,K) 

TSN 0015 DX=XYZ (Nil, 13 - XYZCNL, 1) 

TSN 0016 GO TO 300 

ISN 0017 100 NH=NFC2» J»K) 

ISN 0019 DX=XYZCNH,1) - XYZCH,13 

ISN 0019 GO TO 300 

ISN 0020 . 200 NL=NFC1I-1§ J,K) 

ISN 0021 OX=XYZ (N » 1) - XYZCNL, 1) 

ISN 0022 300 CONTINUE 

ISN 0023 IF C J -EO. 1) GO TO 400 

ISN 0025 IFCJ -EQ. JJ) GO TO 500 

ISN 0027 NH=NFC I* J+1,K) 

ISN 0028 NL=NFCI, J-1,K> 

ISN 0029 DY=XYZCNH,2) - XYZCNL, 2) 

ISN 0030 GO TO 600 

ISN 0031 400 NH=NFC I , 2, K3 

ISN 0032 OY=XYZ (NH, 2) - XYZCN.2) 

ISN 0033 GO TO 600 

ISN 0034 500 NL=NFCI , JJ-l.K) 

TSN 0035 DY=XYZCN,2) - XYZCNL, 2) 

ISN 0036 600 A=DX*DY/4- 

ISN 0037 IFCKP -EQ. 1) GO TO 700 

ISN 0039 IF C MP .EQ. 2) GO. TO 800 

ISN 0041 RETURN 

ISN 0042 700 NU=0 

ISN 0043 IF CK .EO. KK) RETURN 

ISN 0045 NH=NFCI, J,K*1) 

ISN 0046 NU=NH 

ISN 004 7 D ? =X YZ CNH » 3) - XYZCN,3) 

ISN 0048 GO TO 900 

ISN 0049 800 N0=0 

ISN 0050 IF C X -EQ. 1) RETURN 

ISN 0052 NL=NF<I, J,K-1) 

ISN 0053 ND=NL 

ISN 0054 OZ=XYZCN,3) « XYZCNL, 3) 

ISN 0055 900 C=A/DZ 

ISN 0056 RETURN 

ISN 0057 END 


3rj 

s 

i 

m 

e 

3& & 

r* n* 


4<7 


< 
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0S/360 "FORTRAN 11 


0004 

0005 

0006 

0007 

0008 
0009 
0011 

0013 

0014 

0015 

0016 
00 J 7 
0018 

0019 

0020 
0021 
002 2 
0023 
0025 

0027 

0028 

0029 

0030 

0031 

0032 

0033 

0034 

0035 

0036 

0037 
0039 
0041 
00 4 2 
0043 
0045 

' 0046 
004 7 

0048 

0049 

0050 

0052 

0053 

0054 

0055 

0056 

0057 


COMPILER OPTIONS - NAHE= MAIN, OPT =02, LINECNT=60 . S IZE=OQ0OK, 

SOURCE,EnCOIC,NQLIST,NQOECX,LQAQ,HAP,N9EOIT,NOID,NOXREF 
' SUBROUTINE CEUCMP,C,A3 

l COMHON/BK/IIM,IIP,JJM,JJP,KK,NTQT,IV* JV» II ,JJ*H*N»VPC303» 

1 XYZC 2080 , 33 , VVC 30 , 20 , 103 , XPC 303 «XHC 103 ,YPC 203 jYHClQ 3 iZZC 103 * 
2 XXC 4 0 3 » YYC 3 G >, ILX , IUX , KUK , MBC, KB 0 , V RF, NFPS , SKPRFL, SXPLST 
» C 0 H M 0 N / C C C / C N , C S , C E , C W , C U , C D , C C * N N , N S , N E . 1 1 W , If U , N D 

I NFCIX, JX,KX 3 = IX+ II*(JX-U * II* JJ *CXX -13 

» A= 0 . 

r c=o. 

1 CALL FINDCIfJ.KT 

) IFCJ-E0-13 GO TO 100 

IFCJ -EQ- J J 3 GO TO 200 

I NM=NFCI, JH,K3 


- XYZCN.23 


- XYZCNL, 33 


XYZCN,33 


NL=NFCI, J-1.K3 
DY=XYZCNH, 23 - XYZCNL, 23 
GO TO 300 
NH=7iFC I * 2*10 
DY=XYZCNH,23 - XYZCN.23 
GO TO 300 
NL=HFCI » JJ-I,K3 
0Y=XYZCN,23 - XYZCNL , 23 
CONTINUE 

IF CK .EQ. 13 GO TO 400 
IFCK . EO* KK3 GO TO 500 
NH=NfCI, J,K>13 
NL=NFCI, J,X-13 
DZ=XYZ CNN, 33 - XYZCNL, 31 
GO TO 600 
NH=NFCI. J,23 
DZ=XYZCMH,33 - XYZCN,33 
GO TO 600 
NL=NFCI» J» KK-13 
DZ=XYZCN,33 - XYZCNL, 33 
A=DY*DZ/4. 

IFCNP -EQ. 13 GO TO 700 
IFCMP .EQ. 23 GO TO 800 
RETURN 
NE = 0 

IFCT.EQ.II3 RETURN 
N 11 = NFC 1+1, J,K3 
NE-NH 

UX=XYZCNH,13 - XYZCN.13 

GO TO 900 

NM=0 

IF C I .EQ. 13 RETURN 
NL=NFCI-1. J.K3 
NW = NL 

DX = XYZCN ,13 - XYZCNL s 13 

C=A/DX 

RETURN 

END 


-i * 


i 


G S/360 FORTRAN II 


Yl 1. 

il.H 

C JUH 7 4 

j G S/360 FORTRAN II 


COMPILER 

OPTIONS “ NAHE = HAIN^OPT=02»LINECNT=60»SIZE=OOQQK» 

SOURCE » EBCDIC, NOLlSTs NOQECK, LOAD, MAP, NOEOIT.NOIO.NOXREF 

ISN 

0002 


SUBROUTINE CNSCHP»C,A) 

I SN 

0003 


COMHUN/BK/IIM, LIP , JJM, JJP.KK .NTOT.IV, JV,II , JJ,H,N,VPC30) , 
1XYZC2080,3), VVC30,20,10),XPC3G),XMCl6),YPC203,YHa0) ,ZZC10), 
2XXC40) , YYC30) ,ILX, IUX , X.UK , HOC , HBD, VRF, NFPS , SKPRFL, SXPLST 

ISN 

0004 


COMMON/CCC/CN»CS,CE»CW»CU,CO,CC»NN,NS»NE,NW, NU*NO 
NFC IX, JX,KX)=IX*I I*CJX-L 3 *11* JJ*CKX-13 

ISN 

0005 


ISN 

0006 


A = Q . 

ISN 

0007 


C=0. 

ISN 

0008 


CALL FINDCI,J,K) 

ISN 

0009 


IFCI.EO.n GO TO 100 

ISN 

0011 


iFCI.EQ.m GO TO 200 

T SN 

0013 


NII-NFC IH.J.K) 

T SN 

0014 


NL=NrCI-l, J,K) 

ISN 

0015 


OX=XYZCNH, 1) - XYZCNL, 1) 

I SN 

0016 


GO TO 300 

ISN 

0017 

100 

NH=NF C 2 , J, K) 

ISN 

001 8 


UX = XY7.CNH,1) - XYZCN.l) 

I SN 

0019 


GO TO 300 

ISN 

0020 

200 

NL=NFCII“1,J»K) 

I SN 

0021 


OX- XYZCN ,1) - XYZCNL, 1) 

ISN 

0022 

300 

CONTINUE 

I SN 

0023 


IFCK.EQ.l) GO TO 400 

ISN 

0025 


IFCK.EO.KlOr GO TO 500 

ISN 

0027 


NH=MFCI» J.X+13 

ISN 

0028 


NL=NFCI, J,K-1) 

T SN 

0029 


DZ=XYZCNH,3) - XYZCNL, 33 

ISN 

0030 


GO TO 600 

ISN 

0031 

400 

NH=NFCI, J,2) 

ISN 

0032 


0Z=XYZCNH,3) - XYZCN.3) 

ISN 

0033 


GO TO 600 

ISN 

00 34 

500 

NL-NFCI, J,KX-1) 

ISN 

0035 


DZ=XYZCN,3) - XYZCNL, 3) 

ISN 

0036 

600 

rt=0X*DZ/4. 

ISN 

0037 


IFCMP . ED-13 GO TO 700 

ISN 

0039 


IFCMP.EQ-2) GO TO 800 

ISN 

004 1 


RE I URN 

ISN 

0062 

700 

NM = 0 

ISN 

00 4 3 


II'CJ-EQ.JJ) RETURN 

ISN 

004 5 


MH=NF(Is JU.X) 

ISN 

0046 


NN = NH 

ISN 

004 7 


DY = XY ZCNH, 2) - XYZCN * 2) 

ISN 

0048 


GO TO 900 

I SN 

004 9 

800 

NS = 0 

ISN 

0050 


IF(J-EQ-l) RETURN 

ISN 

0052 


NL=NFCI,J-1,K) 

ISN 

0053 


NS=NL 

ISN 

0054 


OY = XYZ CN ,23 - XYZCNL, 2) 

ISN 

0055 

900 

C=A/DY 

I SN 

0056 


RETURN 

ISN 

0057 


END 



TSN 


ISN 

ISN 

ISN 

ISN 


ISN 


ISN 

ISN 

ISN 

ISN 

ISN 

TSN 

ISN 

ISN 

ISN 

ISN 

ISN 

TSN 

ISN 

ISN 

ISN 

ISN 


ISN 
TSN 
ISN 
ISN 
ISN 
ISN 
ISN 
ISN 
ISN 
I SN 
ISN 
ISN 
TSN 
ISN 


COHPI 


0002 


0003 

00 0 4 

0005 

0006 


0007 


0008 

0009 

OOLO 

0011 

0012 

0013 

0014 

0015 

0016 

0017 

0018 
0019 
0021 
0022 

0023 

0024 


0025 
002 6 
0027 
0020 
0029 
0031 
0033 

0035 

0036 

0037 

0039 

0040 
004 I 
0043 


LElR OPTIONS - 


C 

c 

c 

c 

c 


c. 

c 

c 

c 

c 

1001 

2001 


3000 


200 

C 

c 

c 

c 

c 


SUUROUTT 


N'AME = HAIN,0PT=02,LINECNT=60,5TZE=0000K 
SOUR CE, EBCDIC, NOLI ST, NOB ECK.,LUAU , HAP, NOE 
NE FIELD 


NOEDIT, NQID.HOXRSF 


CONSTRUCTION OF COEFFICIENTS CHATRIX ELEMENTS') 

IN LINEAR DIFFERENCE EQUATIONS 
SOLUTION BY OVEP RELAXAT ION 

C0MM0N/3K/IIM, IIP, JJM, J JP , KX . NTOT , IV, JV, II , JJ,M, N, VPC303 , 
1XYZC2 08 0.3). VVC3 0, 20,10), XPC30), X MCI 0) , YP C20) , YHC1 0) , ZZC1 0) , 
2XXC40) , Y YC 30 ) , I LX , IU X , KUK , MBC » MBD . VRF , NFPS ,SKPRFL,SKPLST 
COMMON/FLQ/XC2O8O,2),COEFC208O,7),INDXC2OB0,6),SKPCO 
COMMON/CCC/CN,CS,CE,CW,CU,CD,CC,NN,NS,NE,hW,NU,ND 

t 

INTEGER GQ,ON/*NQRT' / , OS/ * SOUT * / , OE/ * E AST * /, OW/ « WEST * / , 

1 OU/'UP 1 /, OD/* DOWN * / 

ASSUME ASYMPTOTIC MONOPOLE AT INFINITY 

ALPHAF(UUU)=A8SCUUU/RS) 

NDO=POSITI VE FOR DIAGNOSTIC OUTPUT 
ND0 = 0 

WRITECM, 1000) 

FORMAT C1H1/18HCFIEL0 CALCULATION) 

WRITECM, 2000) 

FORMAT C////1X, I7HC0EFFICI ENT ARRAY) 

X0 = .5*XP CIV) 

Y0=.5*YPCJV) 

Z0LD=0 . 

DO 600 N=1 , NTOT 

RS=CXYZCN,1)-X0)**2 +CXYZCN,2)-Y0)**2 +XYZCN,3)**2 
CALL FIND Cl , J, K) 

IFCZZCK) .LE. ZOLD.AND.N.GT.l) GO TO 200 
ZQLD=ZZC K.) 

WRITECM, 3000) K,ZZCK) 

FORMAT C //1X,2HZC,I2,2H)=,F6.3/ 

1 12X,1H0,17X,1HS,17X,1HW,17X,1HC ,17X , IUE , 1 7X , 1HN , 17X , 1HU) 

l» L U« 


MODIFICATION TO SOLVE HELMHOLTZ EQUATION USING LINEARIZED SPACE 
CHARGE- HELM = 0 E~3 Y E-LENGTH-LIKE PARAMETER. CASSUHES POTEN- 
TIALS ARE DIMENSIONLESS) 


HELH=0 -0 
V0LSQ=1. 

DO 300 MP=1 , 2 
CALL CNSCMP.C, AREA) 

IF CMP-EQ.l) OQ=ON 
IF CMP.EQ.2) 00=0S 

IF CNDO-GT-O) WRITE CM.8B8) N, I, J, K , 00 , ARE A,C 
88 8 FORMAT C1X,I8HN,I, J.K.OO, ARE A, C = , 14, 2X , 313 , IX , A5 , 1P2E16 .4) 

IFCC.GT.O.) GO TO 250 
YYY=XYZCN,2)-Y0 
ALP11A = ALPHAFCYYY) 

IF CNOO.GT.O) WRITE CM, 999) N , I , J, X , ALPH A 
99 9 FORMAT C1X, 14HN , I , J , K , A LP HA = , I 4, 2X ,3 13, IP El 6. 4) 


0 0 4 4 
0045 
0047 
004 9 

0050 

0051 

0052 
00 5 3 
0055 
0057 
00 59 
0060 
0062 

0063 

0064 
0066 
0067 
0069 

0071 

0072 

0073 

0074 

0075 
0077 
0079 
00Q 1 
0082 

0084 

0085 
0087 
00 8 0 
0090 

0092 

0093 
009 4 
0095 
0097 
009 8 


CC=CC+AREA*ALPHAFCYYY9 
IFCMP.E0.L9 CN=C 
IFCHP.EO.29 CS=C 
CONTINUE 

VOLSO=VOLSQ*AREA 
DO 400 H. P = 1 » 2 
CALL CEWCHP, C, APEA9 
IF CMP. ED. 19 00=0£ 

IF CMP. ED. 29 00=0 W 
IF CN00.GT.09 WRITE CM, 8889 
CC=CC+C 

IFCC.GT.0.3 GO TO 350 
XXX=XYZCN,19“X0 
ALPHA=ALPHAFCXXX9 
IF CMDO.GI.09 WRITE CM, 9999 
CC=CC*AREA*ALPHAFCXXX9 
IFCMP-EQ.13 CE=C 
IFCMP . EO . 23 C W=C 
CONTINUE 

VOLSQ=VQLSQ*AREA 
00 500 MP=1,2 
CALL CUD CMP , C , AREA9 
IF CMP. 60. 19 DO=OU 
IF CMP. EO. 29 OD=QD 
IF CNDO.GI.03 WRITE CM, 8889 
CC=CC*-C 

IFCC.GT.O..OR.CC.EO.O..ANO. 


CM, 8889 N,I,J,K,QO,AREAfC 


.K, ALPHA 


8889 N,I, J,K,00,AREA,C 


IFCC.GT.0..0R.CC.E0.0..AN0.MP.E0.299G0 TO 
ALPHA— ALP HA F CXY Z CM, 399 
IF CND0.GT.09 WRITE CM, 9999 N , I , J , K, ALPHA 
CL = CC+-AREA*ALPHAFCXY2CN,399 
IFCMP. EQ. 19 CU=C 
IFCMP. EO. 29 CD=C 
CONTINUE 

VOLSO=VQLSQ*AREA 
VOL=SQRT CV0LS09 

IFCHELH.GT.0.9 CC=CC+V0L/HELK**2 

CALL ARRAY 

CONTINUE 

CALL RELAX. 

RETURN 
E M 0 


►* I* * * 


I 


* ->uw 


i ua i a«ri n 


DATE 


COMPILER OPTIONS - NAH£ = H AI N , OP T=0 2 t LIMEC NT=6 0, S IZE=QQ00K s 

SOURCE, EBCDIC, NOLIST, NODECK, LOAD* HAP tNOEDIT, NO 10, NDXREF 
ISN 0002 BLOCK. DATA 

TSN 0003 C0MH0N/CD/PV0LTS tXMACH 

ISN 0004 REAL P AR Tl/ ' IONS ‘ * 

ISN 0005 END 


,DENST.NN« PARTCLC2),PART1C2),PART2C2) 
'/jPARTZ/^ELEC'i’TRON*/ 


Appendix B 


Computer Programs: Electric Fields Produced by 

Cloud-to-Ground Lightning Flashes 

The following four pages contain a listing of the computer programs 
written to compute the electric field produced on the ground as a function of 
time and distance from "ground zero" by the charges associated with a cloud- 
to-ground lightning flash. This program was written by Jerry L. Bohannon. 


T.ITL ClOUD-rC' GROUND SIMULATIIOA 
BATCH 

UE= STROKE 
JBUG 

I MPL I C I T INTEGER +2 ( I-N> 
lOlMENSION RS I (2, 10*-. RSI S, ( 2, 1 0 I 
•OAT A TPIE/5- 5 60 6 2E- U /. TIP HV2.0E-? / 

•DATA IC ARDS/'C"/ 1 , MERM/“M"/*:i Y/" Y~/,L'N/"N"/ 

data ima/X" MU 5-v. icy / x~:i ci6T /. ib el / sro 7 :7 '/„pb g/x" iEio*v • 

DATA rR 0/ X" 1!0 11*"/ ,2-GR / XT1C1 2" /, I V EV X" 101 3"/ , f.«L/X"lD14“/ 
•DATA I A7N/)C" , 0 E" /» IA7.F/X*";0 F* /:. NULL/ X"00" /» IHOMEVrOB-/ 

•DATA IBGV/)T'1E13''/ 

DATA RS I( 1, 1 >/0.0/.RS I<2 ,1>/0.0/', RSI< 2* 10)/ 0-0/ 

DATA PIE/ 3- 14 15 9 2 6** *. 

•DATA RHG/2-0 E-9/. 

CROOT= 1-/3- ' 

WRITEIlA.l) I 
FORMATS "1") 

DO 1000 1=1, 32000 
JCC ,K=I 

WRITE (14,4) IIA7N. WA. IB8L.IA7F 

lF0RWAT<2A2-”LICHrft£X-nG SQ.LT S I MJLA Til GN-+- .R0l’*,.2A2) 

} WRITE! 14,11)’ IBL. [GR 

L FORMAT! A2."READ (DATA FROM CARDS .'OR TE-RM'I N AL",-A 2 ) 

READ! 15. 12) : I WHERE 
! '.FORMAT ( A 1 > 

I :IF< CWHERE.EQ. ICAR OSfSGOTO- 50 

•rFMWHERE.EQ. IT£RM> UGCTO. 70 

WRITE! 14,14) IRO, IRG 
. FORMAT! A 2." TRY AGAIN" ,A2 > 

GOTO 10 

i •READ(1,51,END=999 ) -TO,OGL ,QSa.,VSL ,fVRSJ( CRSHUJ), 1=1,2) ,J=2 
FORMAT < 5< F6. 0.2X )/8(2F10-.C/ T-) 

IGOTO 90 

WR ITE( 14 , 71 ) IMA 

.FORMAT< A 2," E“NTERT FLOATI NG POINT INITIAL [CONDITIONS F6-0") 
'WR1TE!14.75 )i IBL. ICY! 

FORMAT (A 2," f-N IT UAL .-HEIGHT KM", A 2 ) 

READd-5,731 lYO 
WRITE(I4.76) IBL. IC Yl 

FOR MAT ( A2," CL OUO CHARGE , C0Ul'.A2)* 

■READ! 15,73) IQCL t 
'WRITE ( 14 , 72 ) I8L,rCYI 

•FCRMAT( A2 , "STEPPED LEADER CHARGE COOL’.A2) 

READ( 1-5,73) :OSL r. ORIGINAL 

FCRMAT(F6.0) OF POOR 

WRITE< 14,74) IBL.ICYl 

F0RKAT(A2 ."STEPPED LEADER VELOCITY E5 M/S',A2) 

REAO(15«73> IV SL l 
'WRITE (14, 77)1. IBL, ICYI 

FORMAT! A2." RETURN STROKE! VELOCITY E7 M/S".A2) 

RE AD (15 *73) :VRS " 

WRITE( 14. 80) IBL. ICY! 

FORMAT! A2 ." ENTER 8 -TIMES (MS) AND CURRENTS OCAHP) TO DEFINE 
•o- RETURN STROKE 2F1 0 . 0"/:' " * " . A2 ) 

00 32 J=2.9 

RE AD <15. 31) RSI (1 „J>-.RS I ( 2, J ) 

IF(RSI( 1. J)- LT-0. •) .IGCTO. 73 
CONTINUE 


PAGE is 
QUALITY 


THE"!/ 


SI .FORMAT < 2 F 10. 0 ) 

90 VSi.S= VSL 

RSI'C1.10)=RSI(U9) 

‘DO 10C2 J=l,10 
RSIS! 1. J ) =RS I (1. J> 

RSIS!2.J)=RS 1(2. J) 

RSIC.l. J» = RSI ! U JI/10JO- 
PSU2.J) = RSI ! 2*~J ) * 10.10. 

.0.02 .CONTINUE 
:Oa=-QCL 
QSi-=-CSL 
!Y0S=ir0 

RAD=!0.75*ABS!QCU/*RlH0/PJtE>**CRa0r 

DO *1005 1=2. 9 

A=R5I C 2 • I ) 1 

8=RSI ! 2 . I— I I 

C=RSI (2 .1+1 )* 

TF< AL.CT . B-AND.A-.GT . C ) X<T.T = l 
035- IC0NUNU6 
VRSS= VRS 

00 VSL=— VSL * 1- 0 £5 
■ Y0=Y0*1 OOC-O 
!VRS=VRS*1.0E7 

• OTSL=I . OE-4 ’ 

05 WRITE! 14, 110 ) IBL.ICT 

10 FORMAT! A 2." WHAT ‘.IS RA0UTS",A2) 

•READ!15.73» 'D 

WRITE! 13. Ill I IR'D.IB*GY 

11 FORMAT! 2A2/"U."> . 

DO 1 DO 1 1=1. 32000 

001 K = I 

WRITE! 13T 149 > 

149 FORMAT! l-X.'S-I UNITS"//! 

50 <WRI T£U3,151«r YOS.QCL.QSIL ,VSLS. VRSiS*R'SI S'. D 

51 FOR MAT! IX ."HE IGHT= , F7 — IC ," KM" / 1 X . "0— CL QUD = ",F7.1«" C" / 
SlX.*“Q-LEADERI= ',F6.1,' C," /I X. "V- LE ADSR= ".F6.1," E5 M/tS"/lX. 
S "V-RETURN= % F6. 1. ' E7. M/S" /IX." RETURN STROKE MS. KAMP'/ 

Si 2F IQ -4 />/!/ IX." RAOtl US= ",F6-0." M "///)' 

WRITEU3. 152 I IR-G. IRD 

5 2 FORMAT < 1 X. 

S A2.aX."T".15X.'E''.l?X."a;",16X."H'. A2I 
T = 0_0 

slry 2 =i .o/< DF»D*ra*ra') 

SLRY=SQRT!SL RY2 I 

• YC= Y.O-+R AD 

StRQC2=l. v/rD*D-+YC*YC> 

SLRQCL = SQRT!lSlRQC2)-*SLRQ C2 
DI=1.0/D 
X = YO 

EM AX-0.0 
CONTil fiUE 

:slrx 2 =i.c/! d>d*x*x> r 
SLRX=SQRT (SEP X2 I 
SLRX32=SLRX**SLR X2 

E=QSL/TP IE/T.O*! SLRK-SLR Y;)-*SfR CCL * YC/TPI £* ! QCL— QSL*< l.-X/YO* » 
IF! A8S! E ) JST. A8S! EMA'X) ) .E MA X= E 
IF! ASS! E) .LTf.5.0 E4) GOTO 211 
•WRITE! 13.210 ) T, E * X 
C FORMAT! F 16. 7 . F16-0..16X. F16.ll 



215 
21 1 


3C0 


301 


10 


22 

20 


! 1 


2 

0 


1 


FORdiAT* FI 6- 7. F 1'6.0, F16„;5 ,F 16 .1 . 1 1 C *F 16 . 7 , F 16 . 0 ) 
IF< ABS( £) .L 71. 5. C£4. OR.X.IGT.'' . 5E 3 ) D IS 1=1 . VE-3 
T=T~*-DTSL 

IF<X;LT.50. ) GOTO 50 i 

x=ro+vsi *t 

IFCX.LT. J.O I GOTO 50: 

DTS1. = 1 . 0£-4 ' 

: GOTtJ 20 0 
.CONTINUE 
T=T— .OTSL 
WRITE(13*50l ) 

FCR HAT ( IX*' : “) 

ESC=E 

QRC=QSL— QCL t 
PL=~QSL/7D 
SLRT.03= YC*SL RQCti 
•KRNT=i 
R 1=0. 0 
iTR=Q_ 0 
•K01.D= 0 
Q=0. 0 
CONTINUE 

■CAIU CURENT<tRSI.O.»T, TR.IR I, KJR AT .XGLDJ 
IFTRI.LE.UO > GOTO 6 2 j 
T=T-*DT 
Y=WIS*TR 

ifcy.gt.yoj rGaro 522 

P=Q/Y 

:SLRYR=1 _C /sort* cx+D*r,* Y > 

E=ESt,+ P* ( Dl— "SLRY R J/T-P IE 

IF<A8S*E).GT.ABS(EM4X)I ,El'AX=E 

IF< TQT.LT.KRNT..AND. ABS<0) XT.5.CE4 > GOTO 51"’ 

WRI TE( 13* 21 5 > T».E.O^Y.XR.NT ,)TR,RI 
GOTO 510 
'WRITE* 1 3*501 ) 

CONTINUE 

P=Q/.YO 

IF* P -GT.PLJ IGOTtr 572 

E=ESL+P* (Dl—~SLRY I/TPIE 

IF* ABS(E) .GT . A8S< EMAIX )> .Ef'AX=£ 

IF* IQT.LT.Kff.NT. AND. ABS(E') XT .5.CE4 > GOTO 521 
'WRITE Cl 3*215 ) ' T *'E« 0-.IY 0*KIRAT*!TR*RI 
•CACU. CURENHIRSI «T3.DT,TR*|R I.WRAT. KOLD» 

IF(R I.'LE .0.0 ) GOTO 6. J 
T=T+DT 
GOTO 52 C 
■WRITE (13* 501 > 

ICGNTI NUE 
OR S=;Q-*-QSL 

IF*ORS-GT.ORtC> GQTDr.600 

E=ESL + PL**0 I;-SLRsY ) /■ T.P IE*)Q RS *!SLR Y 03 / TP IE i 
IF* A8S< E ) ..GT .A8S* EM A‘X ) ) .EHA.*=E 
IF* IQT. LT.KfTNT. AND. A85* 0) XT.5.CE4 > GOTO 571 
•WRITE* 13,215) T.X-Q.IYC. KIRNT.lTR* RI 
•CALL’ CURE NT (tRSI*O.OT,TR,IR I, K'R AT * XOLD) 

IF< R l-LE .0.0 ) GOTO 6,.C 
T=T-*DT 
GOTO 570 

WRITE* 1 3* 599 ) Q * EM A X 


0/7 POOR <s 


Quality 


» 


?9 F0RHAT(//1X,'"CRT= “ » F 10„i4 »' C-,5X,'£MAXs " , El 2 . 4- " V/HT) 
•WRITE! 14,601) I BG, IB L, ICR. I EEL 
H FOR 'MAT! 2A2»~ DO YOU -WANT IAAOT HER R AT) [U S' • 2 A 2 ) 

•2 REAO( 15. 12) 1 1 AO 

IFflAD.EQ.ir ) GOTO 1 .5 
I F< IA0 . EQ . I W) GOTO 651 
’WRITE! 14, 14) IRO,IGR. 

GOTO 60 2 

>0 WRITE! 14,651-) I B1. » I GR 

>1 FORMAT! A2.-DO YOU WANT AINCTHER EVENT, A2) 

>2 REAO( 15, 12) IE 

IF! TE-EQ. IV >: GOTO 13 
IF! IE. EQ. IN) GOTO 95., 

•WRITE <14.- 14) IRD^IGRi • 

IGQTO 652 

19 UR IT E< 14-, 99 8. ) IRT) 

18 F OR MAT ( A 2',“ NO MORE CAROS") 

10 WRITE(I4,951r) ItTA,IGR 
;i FCRliAT ( A2 •"“BND OF PROGRA M', A 2 ) 

STOP 1 
END 

AB= CURENT 
BUG 

SUBROUT INE CUREfJT! RSI ,Q,!DT* T,R,RI ,.XIRM,K OLD) 

iMPoicrx integer *z ci-n) 

D ItiE NS I ON RSI <2, 11) . 
iF!K0LD. EQ-KRNT it GOTO 50 
TAU=RSKl,KRfNT*L)-RS!I<l,lKfiNr ) 

IFTTAU.LE-O.O) GOTO il 00 
IF<TAU.LE.l- OE-5 ) DT = 0.5*E-6 

IF<TAU.GT.U C Er-5- A NO.. TALC. L£ . 1 ,0E— 4 ) OT=l..'E-5 
I F< TAU. GT .1- 0 E-4 - AND - T AtX . LE .. 1 -0E-3 ) DT=1. ,£-4 
IF<TrAU.GT-l- 0E-3- AND - T AU.. LE - l .OE-2 ) -DT= 1 , :E-3 
IF!TAU.GT.l- l'E-2 , AND .TAIL. LE .'l .OE-1 ) ,DT=1..'E~2 
IF<TAU. GT.l. CE-U DT.= C-0-25 
DEL !R S 1 <2-, KRNT.^*- D-r RSI (72 R N T ) )VTAU 
I F<t)T EQ . C .0 ) RI 1 = 0 - 0 
<TRR=TR+DT 

IF<TRR-GT.RSI(1,KRNT^ 1) ). CT=RSI( 1, KRNT+ 1 >-TR 

RI=RI-*DELI*DT 

!R I2=R I 

Q=Q-*DT * CR I2-+R 11 >/2- 
'TR=TR-*DT 
KCLD=KRNT 

[ F! IR- GE . RSI < l,KRNT-*-l>) |KRNT-=KRNT+ 1 

R 1 1= R I 2 

RETURN 

CONTINUE 

R 1=0-0 

RETURN 

END 

END 



Appendix C 

Computer Output Listing: Cloud-to-Ground Lightning Flash Density 

The following seven pages are the computed output from the program that 
calculates the lightning flash density (cloud-to-ground) from the monthly 
thunderstorm days using the Pierce Conversion. This program, written by Jerry 
L. Bohannon, uses the Normals, Means and Extremes data from "Local Climatolog- 
ical Data — Annual Summaries for 1977" published by the National Oceanic and 
Atmospheric Administration, Environmental Data Service, Asheville, North 
Carolina (available also on magnetic tape). 
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[43 j THE tOTTOH OF THE MEASUREMENT AREA IS AT GRCVHD „EvS» f H-IwS THE TC» 

[33 a TS ‘UP 1 METERS HIGH, 
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£133 FfP+1 

£1?3 H(RfWID)/LOOPL 

£203 Gf-Gfi 

£213 H(OfHI)/LOCFH 

[223 TRYHTRftl 

£233 

[243 'THIS IS RUH HWMSES 1 , ( t TRY), DATE 

£253 'THE CALCULATED VALUES CF THE ELECTRIC fOTPITJa*, IK vOLTS, ARE SHOWN DELCW , • 
r"4i 1 * 

£273 " 

[23] AARf-POTf-gPOT 
[293 5PO flV C2013 * 
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£373 C ? Fp 1 
£333 5f: sv ::oi: 

[39] jj si3H=0//0 
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” J i • ii 

[42] XFOT 

[433 'THIS IS RUH HUMSER • , ( t TRY) , DATE 
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[133 ■){Rl<WID)/LOOFL 
[193 

[2C3 ^OiiHl)/LCC?H 
[213 RROTrS'F.OT 

[223 7RTi.TR?i + i 
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C243 C^KO 
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[13 A THIS FL'IfCTIC'l COMPUTES THE ELEC'SIC POTEHTthu AT S ti POINT, {L,rt) 
[23 A I» THE SPACE ARCLTO THE AROftf O'? FIVE HLLDOASIS, , 

[33 

H3 WK?,33t\H)WA 
[S] «*3e«*9 

[63 AfrU-KJ 
C73 XU* 

[83 uit-N/o 

[93 »lOOPJHK«.({)USxXJ)t2)*KAf(A-LlCZ3)i2 

[103 MK(H*SXiM)*2)*KA 

CllJ Ut[I3(-+/-(WTB2JSO)X»{HMfCH)*C*5 

[123 IM+I 

[133 ■*( * (H-f J J/SBLOOP 

EH3 UHu.n00000*H)H/u: 
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[13 d THIS FUNCTICH COMPUTES THE ELECTRIC POTENTIAL AT AH I P0IH7, (L,H), 

[23 fl DUE TO THE CHARGED PROTECTION WIRE ABOVE THE BILLBOARD, THIS WISE IS 
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[133 oximrt 

[143 L05P}«Hl(.((H-SlX!iJll*2)+HAlH( i l-UIl[I3)t2 

C153 0M1(.((K + SIX},J1)X2)+HA1 

[163 JXl£X3B+/-(i.w,s2;e3iii*sHi*:jii;,o,5 

[173 IhTtI 

[183 ♦tXiHtl)/UOOP 

[153 ffr/Uii 
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12,42 

TN 

OAK RIDGE 

52,81 

12.71 

TX 

ABILENE 

4 1.73 

7.65 

TX 

AMAR I LLO 

4 e , d 1 

12.25 

TX 

AUSTI N 

40,8 1 

6,7 1 

TX 

6R0WNS V I LLE 

2-4,34 

3.72 

TX 

.corpus christi 

30 a ^6 

4.69 

TX 

DALLA 5=>PT WORTH { DF w ) 

45, 12 

b.04 

TX 

DALLAS (LOVE FIELD) 

AC. 15 

0.62 

T X 

DEL RIG 

35.30 

6 .39 

TX 

EL PASO 

3 5 , a 9 

8,95 

TX 

HOUSTON 

0 7.50 

17.92 

TX 

LUdBOCK 

4 5.32 

10.12 

TX 

MI Du AND='ODcSS A 

Jo. 45 

0.76 

T X 

PORT ARTHUR 

o 4 , 1 7 

to .7 1 

TX 

SAN ANGELO 

3 o , bO 

f ,32 

TX 

SAN ANTGNIC 

3o. 30 

5 « 7 4 

TX 

VI CTCR I A 

aa.13 

4.91 


f 
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DAYS DENSITY 

(NO. /YEAR) (MO./YR ./KM S ) 


TX 

WACO 

45,44 

7,52 

TX 

.V I Ch I T A FALLS 

4<3 1 OD 

4. 3C 

TT 

JOHNSTON ISLAND 

4, 07 

1.17 

TT 

KOR-A ISLAND 

3 6 , •> 5 

5.4c 

TT 

KwAjALEIN ISLAND 

5,76 

1 . 7B 

TT 

MAJURO ATOLL 

16. 52 

2.55 

TT 

PDNAPE ISLAND 

2 6 o 0 4 

i .97 

TT 

TRUK ATOLL 

l 5 9 4 <L 

2,62 

TT 

wAKE ISLAND 

o.93 

1.30 

TT 

YAP ISLAND 

1 a* 03 

2 .46 

UT 

MILFORD 

32 o U c 

7.32 

U T 

SALT i_AK£ CITY 

3 5.25 

0.64 

UT 

■wEN DO V E B 

2 9 , 0 C 

5.77 

VT 

BURLINGTON 

94, 94 

4 . o 3 

VA 

LYNCHBURG 

4 0,50 

9.13 

V A 

NORFOLK 

37. 07 

7.49 

V A 

RICHMOND 

3 5* 75 

7.63 

V A 

ROANOKE 

37. 30 

b. 1 2 

w A 

OLYMPIA 

4i| 6 5 

l .24 

vy A 

SEATTLE (APT) 

7. 27 

1.60 

W A 

SEATTLE (CITY) 

o« 06 

1.46 

W A 

SPOKANE 

1 0. 50 

1 .74 

Uf /\ 

ST AMPEOc PASS 

7*29 

1 . 29 

*V A 

WALLA WALLA 

11.26 

1.51 

w A 

YAK IMA 

6.9 0 

1 .25 

W V 

o£CKL£ Y 

45, 71 

10.97 

« V 

CHAPLSST ON 

43.37 

9.42 

w V 

ELKINS 

4 4 » 9 .1 

10,33 

« V 

HUNTINGTON 

44,36 

9.57 

W V 

PARKERSBURG 

4 4 . 3 0 

9,91 

W I 

GREEN dAY 

34. 79 

6,49 

w I 

•LA CROSSE 

40.15 

6 ,29 

w I 

MADISON 

4 0.62 

7.96 

* I 

MILWAUKEE 

35.3 1 

C • 4 o 

« Y 

CASPER 

34.26 

7, 65 

W Y 

CHEYENNE 

4 9,66 

15.41 

«v Y 

LANDER 

3 1.71 

7. 05 

wY 

3HER I DAN 

3 5.59 

9,0 3 


*■■}■ poc*K 




